Difference between revisions of "Physics 212, 2019: Lecture 5"

From Ilya Nemenman: Theoretical Biophysics @ Emory
Jump to: navigation, search
(An example of a model: Solving for an equilibrium position, harder version)
Line 10: Line 10:
 
;Model building: As before, the equilibrium for the system will be given by <math>\sum {F}=0</math>. There are now three Coulomb forces acting on the particle from the fixed charges. Let's say the charges are at positions <math>p_1<x_2<p_3</math>, with charges <math>q_1</math>, <math>q_2</math>, and <math>q_3</math>, respectively. Let's suppose again that the position of the movable particle is <math>x_1<x<x_2</math>. Then the solution of the problem will be given by <math>\frac{q_1}{(x-p_1)^2}-\frac{q_2}{(x-p_2)^2}-\frac{q_3}{(x-p_3)^2}=0</math>. We notice again that the three charges and fixed charge positions are the only things that are needed for the solution. This is equivalent to <math>{q_1}{(x-p_2)^2}(x-x_3)^2-{q_2}{(x-p_1)^2}(x-p_3)^2 -{q_2}{(x-p_2)^2}(x-p_1)^2=0</math>. This is not a quartic equation in <math>x</math>, and solving this problem would be equivalent to finding the root of this quartic equation between <math>p_1</math> and <math>p_2</math>. This is going to be more interesting!
 
;Model building: As before, the equilibrium for the system will be given by <math>\sum {F}=0</math>. There are now three Coulomb forces acting on the particle from the fixed charges. Let's say the charges are at positions <math>p_1<x_2<p_3</math>, with charges <math>q_1</math>, <math>q_2</math>, and <math>q_3</math>, respectively. Let's suppose again that the position of the movable particle is <math>x_1<x<x_2</math>. Then the solution of the problem will be given by <math>\frac{q_1}{(x-p_1)^2}-\frac{q_2}{(x-p_2)^2}-\frac{q_3}{(x-p_3)^2}=0</math>. We notice again that the three charges and fixed charge positions are the only things that are needed for the solution. This is equivalent to <math>{q_1}{(x-p_2)^2}(x-x_3)^2-{q_2}{(x-p_1)^2}(x-p_3)^2 -{q_2}{(x-p_2)^2}(x-p_1)^2=0</math>. This is not a quartic equation in <math>x</math>, and solving this problem would be equivalent to finding the root of this quartic equation between <math>p_1</math> and <math>p_2</math>. This is going to be more interesting!
 
 
 
 
;Model implementation: We need to find roots of a fourth order polynomial. How do we do this? For this we will introduce the Newton-Raphson method for finding a root of an equation <math>F(x)=0</math>. Let's suppose that we have a guess <math>x_0</math> for a root. We evaluate the function, <math>F_0=F(x_0)</math>, and it turns out that the guess is not quite right, <math>F_0\neq 0</math>. We can then write the Taylor expansion for the function <math>F</math> for the values of <math>x</math> around <math>x_0</math>. Namely, <math>F(x)\approx F_0 + G_0(x-x_0)</math>, where <math>G_0</math> is the value of the derivative of <math>F</math> evaluated at <math>x_0</math>, <math>G_0=\left(dF/dx\right)_{x_0}</math>. But then we can get a better estimate of where the root is by solving the equation <math> F_0+G_0(x-x_0)=0</math>. This would give <math>x_1=x_0-F_0/G_0</math> for the next estimate. Iterating this, one will get <math>x_{i+1}=x_i -F_i/G_i</math>. This way we will approach the solution progressively better and better, as we iterate more. At some point we will need to stop the iteration. Usually, this is yet another parameter that one adds to the initialization in the experiment: the tolerance <math>\epsilon</math> for finding the solution. That is, one interrupts the recursion when <math>\left|x_{i+1}-x_i\right|<\epsilon</math>. Note that the algorithm will only find a root if it is initialized with the initial guess in its vicinity. So if one wants to find a root, one needs to have good guesses. And to find all roots, one needs to have initial guesses near all of the roots. Even then the solution may not exist -- there is no guarantee that the the iterations converge; they can enter cycles or diverge instead. '''Your work:''' with the provided code implementing the Newton-Raphson method, try to find '''all''' roots of the following functions <math>F(x)= {\rm atan} 3x</math> and <math>F(x)=\frac{1}{3+x^2}-\frac{1}{6}</math>. Do you need to modify the code to deal with convergence of the iterations for some conditions? Can you understand why some of your attempts converged to a solution, and some did not? '''Don't forget to submit your code.'''
+
;Model implementation: The analysis above suggests that the entire model consists of solving for roots of a fourth order polynomial. How do we do this? For this we will introduce the Newton-Raphson method for finding a root of an equation <math>F(x)=0</math>. This is the first of many computational algorithms that we will introduce in this class, and it allows to find roots (zero crossings) of more-or-less arbitrary functions, not just of polynomials.
 +
 
 +
:Let's suppose that we have a guess <math>x_0</math> for a root of a function <math>F(x)</math>. We evaluate the function, <math>F_0=F(x_0)</math>, and it turns out that the guess is not quite right, <math>F_0\neq 0</math>. We can then write the Taylor expansion for the function <math>F</math> for the values of <math>x</math> around <math>x_0</math>. Namely, <math>F(x)\approx F_0 + G_0(x-x_0)</math>, where <math>G_0</math> is the value of the derivative of <math>F</math> evaluated at <math>x_0</math>, <math>G_0=\left(dF/dx\right)_{x_0}</math>. But then we can get a better estimate of where the root is by solving the equation <math> F_0+G_0(x-x_0)=0</math>. This would give <math>x_1=x_0-F_0/G_0</math> for the next estimate. Iterating this, one will get <math>x_{i+1}=x_i -F_i/G_i</math>. This way we will approach the solution progressively better and better, as we iterate more. At some point we will need to stop the iteration. Usually, this is yet another parameter that one adds to the initialization in the experiment: the tolerance <math>\epsilon</math> for finding the solution. That is, one interrupts the recursion when <math>\left|x_{i+1}-x_i\right|<\epsilon</math>. Note that the algorithm will only find a root if it is initialized with the initial guess in its vicinity. So if one wants to find a root, one needs to have good guesses. And to find all roots, one needs to have initial guesses near all of the roots. Even then the solution may not exist -- there is no guarantee that the the iterations converge; they can enter cycles or diverge instead. '''Your work:''' with the provided code implementing the Newton-Raphson method, try to find '''all''' roots of the following functions <math>F(x)= {\rm atan} 3x</math> and <math>F(x)=\frac{1}{3+x^2}-\frac{1}{6}</math>. Do you need to modify the code to deal with convergence of the iterations for some conditions? Can you understand why some of your attempts converged to a solution, and some did not? '''Don't forget to submit your code.'''
 
   
 
   
 
 
Line 26: Line 28:
 
 
 
;Discussion:
 
;Discussion:
 +
 
==New Python Concepts==
 
==New Python Concepts==
 
*Loops of different types (for and while).
 
*Loops of different types (for and while).

Revision as of 14:02, 3 February 2019

Emory Logo

Back to the main Teaching page.

Back to Physics 212, 2019: Computational Modeling.

An example of a model: Solving for an equilibrium position, harder version

Problem formulation

A charged particle free to move on a rod is put in the between the first and the second of three immovable particles with an electric charge of the same sign affixed to the rod at various points. Where will the charge come to rest?

Analysis
The same analysis as in the simpler version applies here: it's a static, deterministic, continuous space model. We again need the friction for the solution to exist, and the moving particle with come to rest, and the solution will exist.
Model building
As before, the equilibrium for the system will be given by . There are now three Coulomb forces acting on the particle from the fixed charges. Let's say the charges are at positions , with charges , , and , respectively. Let's suppose again that the position of the movable particle is . Then the solution of the problem will be given by . We notice again that the three charges and fixed charge positions are the only things that are needed for the solution. This is equivalent to . This is not a quartic equation in , and solving this problem would be equivalent to finding the root of this quartic equation between and . This is going to be more interesting!
Model implementation
The analysis above suggests that the entire model consists of solving for roots of a fourth order polynomial. How do we do this? For this we will introduce the Newton-Raphson method for finding a root of an equation . This is the first of many computational algorithms that we will introduce in this class, and it allows to find roots (zero crossings) of more-or-less arbitrary functions, not just of polynomials.
Let's suppose that we have a guess for a root of a function . We evaluate the function, , and it turns out that the guess is not quite right, . We can then write the Taylor expansion for the function for the values of around . Namely, , where is the value of the derivative of evaluated at , . But then we can get a better estimate of where the root is by solving the equation . This would give for the next estimate. Iterating this, one will get . This way we will approach the solution progressively better and better, as we iterate more. At some point we will need to stop the iteration. Usually, this is yet another parameter that one adds to the initialization in the experiment: the tolerance for finding the solution. That is, one interrupts the recursion when . Note that the algorithm will only find a root if it is initialized with the initial guess in its vicinity. So if one wants to find a root, one needs to have good guesses. And to find all roots, one needs to have initial guesses near all of the roots. Even then the solution may not exist -- there is no guarantee that the the iterations converge; they can enter cycles or diverge instead. Your work: with the provided code implementing the Newton-Raphson method, try to find all roots of the following functions and . Do you need to modify the code to deal with convergence of the iterations for some conditions? Can you understand why some of your attempts converged to a solution, and some did not? Don't forget to submit your code.


More to be added.


Model verification

Discussion

New Python Concepts

  • Loops of different types (for and while).
  • We started talking about scopes: what happens to the variable i after the for loops finishes executing?
  • Variables vs objects, yet again.
  • Vector math is always faster than element-by-element math -- we showed this in a simple script in class.