Lab 09.sagews

Bisection Search

Number Guessing Game

• The host have a number in mind

• Audiance try to guess the number

• When the auduance guess a number, the host will answer: "too high" or "too low"
  
  
  If this number ranges from 1 to 100, can all of the students in our class participate?

No, if all the students need to participate, the number need to be larger than 2^19

Describe what has happened!

Pseudocode

higherBound = maxRange
lowerBound = minRange
while(Not found)
    1. guess a number within higherBound and lowerBound
    2. if(guess number > target)
           higherBound = guess
       else
           lowerBound = guess
    3. if higherBound - lowerBound < some small range
       exit while loop

Bifurcation

Quantitative Change Induces Qualitive Change

Goal: Find the Point Where Number of Eq. Points Change!

r =0.1
fig1 = plot(0.5*x^2/(1+x^2),(x,0,6),color='red')
fig2 = plot(r*x,(x,0,6),color='blue')
fig1+fig2

How many equilibrium points we have with current r?

assume(x>=0)
curFun = 0.5*x^2/(1+x^2)-r*x
solve(curFun==0,x)
curFun.roots()
len(solve(curFun==0,x))
len(curFun.roots())

Error in lines 2-2 Traceback (most recent call last): File "/cocalc/lib/python2.7/site-packages/smc_sagews/sage_server.py", line 1188, in execute flags=compile_flags) in namespace, locals File "", line 1, in <module> TypeError: unsupported operand type(s) for *: 'function' and 'sage.symbolic.expression.Expression'

How to get the value?

intPoints = curFun.roots()
for i in range(len(intPoints)):
    float(intPoints[i][0])

0.20871215252208009 4.7912878474779195 0.0

How to Search for Bifurcation

For searching first bifurcation (2 equilibrium points to 3 equilibrium points)

• Find r that has 3 equilibrium points, set this r as higher bound

• Find r that has 2 equilibrium points, set this r as lower bound

• guess a new r in this range

• if has 3 equilibrium points => update higher bound

• if has 2 equilibrium points => update lower bound

• if higher bound is close to lower bound => end

How to guess efficiently?

higherBound = 1
lowerBound = 0.1
assume(x,'real')
assume(x>=0)
rList=[]
while(higherBound-lowerBound>0.0001):
    # ha ha do your own lab


float(rGuess)
fig1 = plot(0.5*x^2/(1+x^2),(x,0,5))
fig2 = plot(rGuess*x,(x,0,5))
fig1+fig2

0.25001831054687507

P.S. solve does not always work!

(however, we will not see it in this lab)

r = 0.3
k = 25
solve(r*(1-x/k)-x/(1+x^2),x)
find_root(r*(1-x/k)-x/(1+x^2),0,2)
find_root(r*(1-x/k)-x/(1+x^2),3,5)
find_root(r*(1-x/k)-x/(1+x^2),20,25)

[x == -1/2*(1/27*I*sqrt(9553223) + 6475/27)^(1/3)*(I*sqrt(3) + 1) - 62/3*(-I*sqrt(3) + 1)/(1/27*I*sqrt(9553223) + 6475/27)^(1/3) + 25/3, x == -1/2*(1/27*I*sqrt(9553223) + 6475/27)^(1/3)*(-I*sqrt(3) + 1) - 62/3*(I*sqrt(3) + 1)/(1/27*I*sqrt(9553223) + 6475/27)^(1/3) + 25/3, x == (1/27*I*sqrt(9553223) + 6475/27)^(1/3) + 124/3/(1/27*I*sqrt(9553223) + 6475/27)^(1/3) + 25/3] 0.32789714277901966 3.621997059419165 21.050105797801766

The End!!