1 files changed, 143 insertions, 0 deletions
diff --git a/assignments/hw02/euler_mod.jl b/assignments/hw02/euler_mod.jl
new file mode 100644
index 0000000..65446d1
--- /dev/null
+++ b/assignments/hw02/euler_mod.jl
@@ -0,0 +1,143 @@
+
+using Winston
+
+function doPlot(V_in, k_l, k_p, K_m)
+    
+    dt    = 0.2
+    tlast = 1000
+    
+    iterations = integer(round(tlast/dt))
+    ATPall = zeros(Float64, (iterations, 1))
+    Gall = zeros(Float64, (iterations, 1))
+
+    ATP = 4
+    G = 3
+    for i = 1:iterations
+        ATPall[i] = ATP
+        Gall[i] = G
+
+        dATPdt = 2 * k_l * G * ATP - (k_p * ATP) / (ATP + K_m)
+        dGdt = V_in - k_l * G * ATP
+
+        ATP = ATP + dATPdt * dt
+        G = G + dGdt * dt
+    end
+
+    tall = dt*(0:iterations-1)'
+    figure()
+    hold(true)
+    plot(tall, ATPall, color="red")
+    plot(tall, Gall, color="blue")
+    hold(false)
+    xlabel("Time")
+
+    println("max(G) = ", maximum(Gall))
+    println("max(ATP) = ", maximum(ATPall))
+    # max(G) = 20.97673324771467
+    # max(ATP) = 18.65978074314273
+    title("")
+end
+
+doPlot(0.36, 0.02, 6, 12.0)
+
+function doXY(V_in, k_l, k_p, K_m)
+    
+    dt    = 0.2
+    tlast = 1000
+    
+    iterations = integer(round(tlast/dt))
+    ATPall = zeros(Float64, (iterations, 1))
+    Gall = zeros(Float64, (iterations, 1))
+
+    ATP = 4
+    G = 3
+    for i = 1:iterations
+        ATPall[i] = ATP
+        Gall[i] = G
+
+        dATPdt = 2 * k_l * G * ATP - (k_p * ATP) / (ATP + K_m)
+        dGdt = V_in - k_l * G * ATP
+
+        ATP = ATP + dATPdt * dt
+        G = G + dGdt * dt
+    end
+
+    tall = dt*(0:iterations-1)'
+    figure()
+    hold(true)
+    plot(Gall, ATPall, color="red")
+    hold(false)
+    xlabel("G")
+    ylabel("ATP")
+end
+
+function testStability(V_in, k_l, k_p, K_m)
+    
+    dt    = 0.05
+    tlast = 2000
+    
+    iterations = integer(round(tlast/dt))
+    ATPall = zeros(Float64, (iterations, 1))
+    Gall = zeros(Float64, (iterations, 1))
+
+    ATP = 4
+    G = 3
+    for i = 1:iterations
+        ATPall[i] = ATP
+        Gall[i] = G
+
+        dATPdt = 2 * k_l * G * ATP - (k_p * ATP) / (ATP + K_m)
+        dGdt = V_in - k_l * G * ATP
+
+        ATP = ATP + dATPdt * dt
+        G = G + dGdt * dt
+    end
+
+    #deltaATP = maximum(ATPall[length(ATPall)/2:end]) - minimum(ATPall[length(ATPall)/2:end])
+    #deltaG = maximum(Gall[length(Gall)/2:end]) - minimum(Gall[length(Gall)/2:end])
+    #return (deltaATP, deltaG)
+    return (maximum(ATPall[length(ATPall)/2:end]),
+            minimum(ATPall[length(ATPall)/2:end]),
+            maximum(Gall[length(Gall)/2:end]),
+            minimum(Gall[length(Gall)/2:end]), )
+end
+
+function doBistable()
+    x = 0.1:0.05:1.6
+    varATP = zeros(length(x))
+    varG = zeros(length(x))
+    ATP_h = zeros(length(x))
+    ATP_l = zeros(length(x))
+    G_h = zeros(length(x))
+    G_l = zeros(length(x))
+    for i=1:length(x)
+        (a_h, a_l, g_h, g_l) = testStability(x[i], 0.02, 6, 12.0)
+        varATP[i] = a_h - a_l
+        varG[i] = g_h - g_l
+        ATP_h[i] = a_h
+        ATP_l[i] = a_l
+        G_h[i] = g_h
+        G_l[i] = g_l
+        println((a_h, g_h))
+    end
+    figure()
+    hold(true)
+    plot(x, varG, color="blue")
+    plot(x, varATP, color="red")
+    xlabel("V_in")
+    ylabel("Stability")
+    hold(false)
+    title("Stability plot")
+
+    figure()
+    hold(true)
+    plot(x, G_h, color="blue")
+    plot(x, G_l, color="blue")
+    plot(x, ATP_h, color="red")
+    plot(x, ATP_l, color="red")
+    xlabel("V_in")
+    ylabel("Value")
+    hold(false)
+    title("Values plot")
+end
+