<jupyter_start><jupyter_text>Question 1<jupyter_code># Sample 100*100 points in the complex square
initial_points = []
zi = []
for i in range(100):
    for k in range(100):
        zi.append([0])
        initial_points.append(i/25.-2+(k/25.-2)*1j)
def iterate(initial_points, zi):
    # Function that takes a point list and the zi list
    # and iterates all points using the equation z_i+1 = z_i^2+c
    
    for k in range(len(initial_points)):
        zi_next = (zi[k][-1])*(zi[k][-1])+initial_points[k]
        zi[k].append(zi_next)
    
    return zi
# Iterate over 20 steps
for i in range(20):
    iterate(initial_points, zi)
pl.figure()
for i in range(1000):
    if abs(zi[i*10][-1]) < 100:
        color = 'blue'
    else:
        color = 'red'
    pl.scatter(np.real(zi[i*10]), np.imag(zi[i*10]), color = color)
pl.xlim(-2,2)
pl.ylim(-2,2)
pl.savefig('question1_1.png')
pl.show()
# we remove any non-finite value of the sequences
for i in range(len(zi)):
    while np.isfinite(zi[i][-1]) == False:
        zi[i].pop()
# Now we assign a color to the index of iteration
pl.figure()
for i in range(1000):
    if abs(zi[i*10][-1]) < 100:
        color = 'blue'
    else:
        color = np.arange(0,len(zi[10*i]))
    pl.scatter(np.real(zi[i*10]), np.imag(zi[i*10]), c = color)
pl.xlim(-10,10)
pl.ylim(-10,10)
pl.colorbar(label="Colour")
pl.savefig('question1_1.png')
pl.show()<jupyter_output><empty_output><jupyter_text>Let's do this again for the box -0.1<x<0.1 and -0.1<y<0.1<jupyter_code>50/0.2

# Sample n points in the complex square
initial_points_2 = []
zi_2 = []
for i in range(50):
    for k in range(50):
        zi_2.append([0])
        initial_points_2.append(i/250.-0.1+(k/250.-0.1)*1j)
# Iterate over 20 steps
for i in range(20):
    iterate(initial_points_2, zi_2)
pl.figure()
for i in range(500):
    if abs(zi_2[i*5][-1]) < 100:
        color = 'blue'
    else:
        color = 'red'
    pl.scatter(np.real(zi_2[i*5]), np.imag(zi_2[i*5]), color = color)
pl.xlim(-0.1,0.1)
pl.ylim(-0.1,0.1)
pl.savefig('question1_2.png')
pl.show()<jupyter_output><empty_output><jupyter_text>They all converge! Wow!  Question 2<jupyter_code># We need scipy for this
from scipy.integrate import odeint
def derivatives(vector, t, beta, gamma, N):
    S, I, R = vector
    dSdt = -beta*S*I/N
    dIdt = beta*S*I/N -gamma*I
    dRdt = gamma*I
    return dSdt, dIdt, dRdt
def ODE_solver(gamma, beta):
    N=1000
    init_vector = 999,1,0
    time = np.linspace(0,200,200)
    final_results = odeint(derivatives, init_vector, time, args=(beta, gamma, N))
    return final_results
result1 = ODE_solver(0.1,0.2).T
result2 = ODE_solver(0.2,0.2).T
result3 = ODE_solver(0.1,0.8).T
time = np.linspace(0,200,200)
pl.figure()
pl.title(r'$\gamma = 0.2, \ \beta = 0.1$')
pl.plot(time, result1[0], color = 'blue')
pl.plot(time, result1[1], color = 'red')
pl.plot(time, result1[2], color = 'green')
pl.legend(['Suspected', 'Infected', 'Recovered'])
pl.savefig('question2_1')
pl.show()
pl.figure()
pl.title(r'$\gamma = 0.2, \ \beta = 0.2$')
pl.plot(time, result2[0], color = 'blue')
pl.plot(time, result2[1], color = 'red')
pl.plot(time, result2[2], color = 'green')
pl.legend(['Suspected', 'Infected', 'Recovered'])
pl.savefig('question2_2')
pl.show()
pl.figure()
pl.title(r'$\gamma = 0.1, \ \beta = 0.8$')
pl.plot(time, result3[0], color = 'blue')
pl.plot(time, result3[1], color = 'red')
pl.plot(time, result3[2], color = 'green')
pl.legend(['Suspected', 'Infected', 'Recovered'])
pl.savefig('question2_1')
pl.show()<jupyter_output><empty_output><jupyter_text>BonusWe want the rate of death to be proportionnal to the number of infected and we want to subtract the deaths to the infected count. Something like this could work:$$\frac{dD}{dt} = \alpha I$$$$\frac{dI}{dt} = \frac{\beta SI}{N}-(\gamma+\alpha)I$$<jupyter_code>def derivatives_bonus(vector, t, beta, gamma, alpha, N):
    S, I, R, D = vector
    dSdt = -beta*S*I/N
    dIdt = beta*S*I/N - (gamma+alpha)*I
    dRdt = gamma*I
    dDdt = alpha*I
    return dSdt, dIdt, dRdt, dDdt
def ODE_solver_bonus(gamma, beta, alpha):
    N=1000
    init_vector = 999,1,0,0
    time = np.linspace(0,200,200)
    final_results = odeint(derivatives_bonus, init_vector, time, args=(beta, gamma, alpha, N))
    return final_results
result1_bonus = ODE_solver_bonus(0.1, 0.4, 0.05).T
pl.figure()
pl.title(r'$\gamma = 0.1, \ \beta = 0.4 \ \alpha = 0.05$')
pl.plot(time, result1_bonus[0], color = 'blue')
pl.plot(time, result1_bonus[1], color = 'red')
pl.plot(time, result1_bonus[2], color = 'green')
pl.plot(time, result1_bonus[3], color = 'black')
pl.legend(['Suspected', 'Infected', 'Recovered', 'Dead'])
pl.show()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.18399657680787335, 'fraction_numerical': 0.05412922550278134, 'mean_word_length': 3.9681190223166842, 'pattern_counts': {'":': 0, '<': 24, '<?xml version=': 0, '>': 17, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '4', 'dir': 'jupyter-structured-clean-dedup', 'id': '445442', 'license': 'MIT', 'n_tokens_mistral': 2060, 'n_tokens_neox': 1931, 'n_words': 425, 'path': 'Assignment2/Assignment_2.ipynb', 'repo_name': 'HenriLamarre/CTA200', 'size': '4674'}