CoCalc Public FilesLab 2 / Lab2-turnin.ipynb
Author: Ryan Kien
Views : 64
Compute Environment: Ubuntu 20.04 (Default)
In [ ]:
# Name: Ryan Kien
# I worked on this code with: Fujia Guo, Nushrat Esha, Xavier Herrera

# Please do all of your work for this week's lab in this worksheet. If
# you wish to create other worksheets for scratch work, you can, but
# this is the one that will be graded. You do not need to do anything
# to turn in your lab. It will be collected by your TA at the beginning
# of (or right before) next week’s lab.

# Be sure to clearly label which question you are answering as you go and to


In :
#1
def dde_solve(dde, statevar, delayedvars, history, tmax, timestep):
# Check validity of delays.
if min(delayedvars.values()) < 0:
raise ValueError("This function will not work with negative delays. "
"Consider consulting a fortune teller instead.")
if any(val<timestep for val in delayedvars.values()):
raise ValueError("Time step should not be larger than delay.")

# Set up variables and delays.

delayedvars = delayedvars.items()
dde = dde.subs({v: statevar for v, delay in delayedvars if delay == 0})
delayedvars = [(v, delay) for v, delay in delayedvars if delay != 0]
allvars = [str(statevar)] + [str(v) for v, delay in delayedvars]
delays = [delay for v, delay in delayedvars]

# Set up fast functions.
dde_func = fast_float(dde, *allvars)
history_func = fast_float(history, "t")

# Adjust the timestep if necessary
mindelay = min(delays) if delays else timestep
timestepcorrectionfactor = ceil(timestep / mindelay)
timestep /= timestepcorrectionfactor

# A function to perform history lookups.
def lookup(t):
"""Does a history lookup at each delay from t, stores result in allvars[1:]"""
for i, delay in enumerate(delays):
if t - delay <= 0:
allvars[i+1] = history_func(t - delay)
else:
r = (t - delay) / timestep
n = floor(r)
r -= n
allvars[i+1] = result[n]*(1 - r) + result[n + 1]*r

# Set up for the first iteration.
result = [history_func(0)]
lookup(0)
for t in sxrange(0, tmax - timestep, timestep):
# Compute k1. Note history lookup has already been done.
allvars = result[-1]
k1 = dde_func(*allvars)
# Compute k2.
lookup(t + timestep/2)
allvars += timestep/2 * k1
k2 = dde_func(*allvars)
# Compute k3. Note history lookup has already been done.
allvars = result[-1] + timestep/2 * k2
k3 = dde_func(*allvars)
# Compute k4.
lookup(t + timestep)
allvars = result[-1] + timestep * k3
k4 = dde_func(*allvars)
# Finally, compute the RK4 weighted average.
result.append(result[-1] + (k1 + 2*k2 + 2*k3 + k4)/6 * timestep)
return result[::timestepcorrectionfactor]

In :
#2

2
In :
L = 6
Vmax = 80
a = 0.2
n = 5
h = 400
tau = 0.2

In :
#3
var("x","xtau")
time = srange(0,100,0.1)
delay={xtau: tau}
xprime = L-((Vmax * xtau^n)/(h + xtau^n))*(a*x)
sol=dde_solve(xprime, x, delay, history=0, tmax=100, timestep=0.1)
list_plot(list(zip(time, sol)), plotjoined=True) In :
#4
output=[]
for x in sol:
output.append((Vmax*x^n)/(h+x^n)*x)

In :
#5
list_plot(list(zip(time, output)), plotjoined=True)
#oscillations occur at a very high rate, and stabilize after a few minutes In :
#6a
tau = .11
var("x","xtau")
time = srange(0,100,0.1)
delay={xtau: tau}
xprime = L-((Vmax * xtau^n)/(h + xtau^n))*(a*x)
sol=dde_solve(xprime, x, delay, history=0, tmax=100, timestep=0.1)
list_plot(list(zip(time, sol)), plotjoined=True) In :
#6b
tau = 1
var("x","xtau")
time = srange(0,100,0.1)
delay={xtau: tau}
xprime = L-((Vmax * xtau^n)/(h + xtau^n))*(a*x)
sol=dde_solve(xprime, x, delay, history=0, tmax=100, timestep=0.1)
list_plot(list(zip(time, sol)), plotjoined=True) In :
#6c
tau = .5
var("x","xtau")
time = srange(0,100,0.1)
delay={xtau: tau}
xprime = L-((Vmax * xtau^n)/(h + xtau^n))*(a*x)
sol=dde_solve(xprime, x, delay, history=0, tmax=100, timestep=0.1)
list_plot(list(zip(time, sol)), plotjoined=True) In :
#7
@interact
def equation(newtau=(0.11,1,.01)):
L = 6
Vmax = 80
a = 0.2
n = 5
h = 400
var("x","xtau")
time = srange(0,100,0.1)
delay={xtau: newtau}
xprime = L-((Vmax * xtau^n)/(h + xtau^n))*(a*x)
sol=dde_solve(xprime, x, delay, history=0, tmax=100, timestep=0.1)
return(list_plot(list(zip(time, sol)), plotjoined=True))
#the minimum delay for persisten osciallations is when tau = .20

In :
#8
@interact
def equation(steepness=(1,25)):
L = 6
Vmax = 80
a = 0.2
n = steepness
h = 400
tau =.2
var("x","xtau")
time = srange(0,100,0.1)
delay={xtau: tau}
xprime = L-((Vmax * xtau^n)/(h + xtau^n))*(a*x)
sol=dde_solve(xprime, x, delay, history=0, tmax=100, timestep=0.1)
return(list_plot(list(zip(time, sol)), plotjoined=True))
#n controls the steepness

In :
#9
@interact
def equation(newtau=(0.11,1,.01), steepness=(1,25)):
L = 6
Vmax = 80
a = 0.2
n = steepness
h = 400
var("x","xtau")
time = srange(0,100,0.1)
delay={xtau: newtau}
xprime = L-((Vmax * xtau^n)/(h + xtau^n))*(a*x)
sol=dde_solve(xprime, x, delay, history=0, tmax=100, timestep=0.1)
return(list_plot(list(zip(time, sol)), plotjoined=True))

In :
#10
@interact
def equation(newtau=(0.11,1,.01), steepness=(1,25)):
L = 6
Vmax = 80
a = 0.2
n = steepness
h = 400
var("x","xtau")
time = srange(0,100,0.1)
delay={xtau: newtau}
xprime = L-((Vmax * xtau^n)/(h + xtau^n))*(a*x)
sol=dde_solve(xprime, x, delay, history=0, tmax=100, timestep=0.1)
return(list_plot(list(zip(time, sol)), plotjoined=True))
#5 is the minimum value for where persistent oscillations occur

In :
#11a
@interact
def equation(newtau=(0.11,1,.01), steepness=(1,25)):
L = 6
Vmax = 80
a = 0.2
n = steepness
h = 400
var("x","xtau")
time = srange(0,100,0.1)
delay={xtau: newtau}
xprime = L-((Vmax * xtau^n)/(h + xtau^n))*(a*x)
sol=dde_solve(xprime, x, delay, history=0, tmax=100, timestep=0.1)
return(list_plot(list(zip(time, sol)), plotjoined=True))
#when tau = .3, the minimum steepness value for oscillations is 4

In :
#11b
@interact
def equation(newtau=(0.11,1,.01), steepness=(1,25)):
L = 6
Vmax = 80
a = 0.2
n = steepness
h = 400
var("x","xtau")
time = srange(0,100,0.1)
delay={xtau: newtau}
xprime = L-((Vmax * xtau^n)/(h + xtau^n))*(a*x)
sol=dde_solve(xprime, x, delay, history=0, tmax=100, timestep=0.1)
return(list_plot(list(zip(time, sol)), plotjoined=True))
#when tau = .5, the minimum steepness value for persistent oscillations is 3

In :
#11c
@interact
def equation(newtau=(0.11,1,.01), steepness=(1,25)):
L = 6
Vmax = 80
a = 0.2
n = steepness
h = 400
var("x","xtau")
time = srange(0,100,0.1)
delay={xtau: newtau}
xprime = L-((Vmax * xtau^n)/(h + xtau^n))*(a*x)
sol=dde_solve(xprime, x, delay, history=0, tmax=100, timestep=0.1)
return(list_plot(list(zip(time, sol)), plotjoined=True))
#when tau = .7, the minimum steepness value for persistent oscillations is 3

In :
#11d
@interact
def equation(newtau=(0.11,1,.01), steepness=(1,25)):
L = 6
Vmax = 80
a = 0.2
n = steepness
h = 400
var("x","xtau")
time = srange(0,100,0.1)
delay={xtau: newtau}
xprime = L-((Vmax * xtau^n)/(h + xtau^n))*(a*x)
sol=dde_solve(xprime, x, delay, history=0, tmax=100, timestep=0.1)
return(list_plot(list(zip(time, sol)), plotjoined=True))
#when tau = 1, the minimum steepness value for persistent oscillations is 3

In :
#11e
@interact
def equation(newtau=(0.11,1,.01), steepness=(1,25)):
L = 6
Vmax = 80
a = 0.2
n = steepness
h = 400
var("x","xtau")
time = srange(0,100,0.1)
delay={xtau: newtau}
xprime = L-((Vmax * xtau^n)/(h + xtau^n))*(a*x)
sol=dde_solve(xprime, x, delay, history=0, tmax=100, timestep=0.1)
return(list_plot(list(zip(time, sol)), plotjoined=True))
#when tau = .11, the minimum steepness value for persistent oscillations is 9

In [ ]:
{"ft_lang_label":"__label__en","ft_lang_prob":0.5960523,"math_prob":0.9997967,"size":7394,"snap":"2021-04-2021-17","text_gpt3_token_len":2550,"char_repetition_ratio":0.14763193,"word_repetition_ratio":0.44509232,"special_character_ratio":0.37489855,"punctuation_ratio":0.19814815,"nsfw_num_words":0,"has_unicode_error":false,"math_prob_llama3":0.99990094,"pos_list":[0,1,2,3,4,5,6,7,8,9,10,11,12],"im_url_duplicate_count":[null,null,null,null,null,null,null,null,null,null,null,null,null],"WARC_HEADER":"{\"WARC-Type\":\"response\",\"WARC-Date\":\"2021-01-21T08:46:18Z\",\"WARC-Record-ID\":\"<urn:uuid:acb6eb5e-6a27-4728-8d34-3047a24438b3>\",\"Content-Length\":\"331501\",\"Content-Type\":\"application/http; msgtype=response\",\"WARC-Warcinfo-ID\":\"<urn:uuid:7a73b624-fede-41db-9ee5-cf3cf13eebaf>\",\"WARC-Concurrent-To\":\"<urn:uuid:de42b5d3-f000-4040-a742-927c9688b2df>\",\"WARC-IP-Address\":\"104.22.1.102\",\"WARC-Target-URI\":\"https://share.cocalc.com/share/ea51c556d2a99093c9328195af52c52493f56ede/Lab%202/Lab2-turnin.ipynb?viewer=share\",\"WARC-Payload-Digest\":\"sha1:CRW4P3Y4MMECNAVXUNDPOVRHJPWL6ZOZ\",\"WARC-Block-Digest\":\"sha1:PHH62NZ4VDJCUTG4O4QWJFINMAR5AHXU\",\"WARC-Identified-Payload-Type\":\"text/html\",\"warc_filename\":\"/cc_download/warc_2021/CC-MAIN-2021-04/CC-MAIN-2021-04_segments_1610703524270.28_warc_CC-MAIN-20210121070324-20210121100324-00620.warc.gz\"}"}