import numpy as np
import matplotlib.pyplot as plt

r_A = 0.9
r_B = 0.4
r_C = 0.7
K_A = 5.5e5
K_B = 5.5e4
K_C = 1e5

n_hours = 24
initial_population = 1000
initial_population_B = 10000
pop_A = np.zeros(n_hours + 1)
pop_B = np.zeros(n_hours + 1)
pop_C = np.zeros(n_hours + 1)

pop_A[0] = initial_population
for i in range(n_hours):
    pop_A[i + 1] = pop_A[i] + r_A * (1 - (pop_A[i] / K_A)) * pop_A[i]
    
pop_B[0] = initial_population_B
for i in range(n_hours):
    pop_B[i + 1] = pop_B[i] + r_B * (1 - (pop_B[i] / K_B)) * pop_B[i]
    
pop_C[0] = initial_population
for i in range(n_hours):
    pop_C[i + 1] = pop_C[i] + r_C * (1 - (pop_C[i] / K_C)) * pop_C[i]
    
data_A = np.loadtxt("data_exp_A.txt")
data_B = np.loadtxt("data_exp_B.txt")
data_C = np.loadtxt("data_exp_C.txt")

plt.figure(figsize=(6,3))
plt.plot(pop_A, label="model A")
plt.plot(data_A, label="experiment A")

plt.xlabel("time (hours)")
plt.ylabel("population")
plt.title("Species A")
plt.legend()

plt.figure(figsize=(6,3))
plt.plot(pop_B, label="model B")
plt.plot(data_B, label="experiment B")

plt.xlabel("time (hours)")
plt.ylabel("population")
plt.title("Species B")
plt.legend()

plt.figure(figsize=(6,3))
plt.plot(pop_C, label="model C")
plt.plot(data_C, label="experiment C")

plt.xlabel("time (hours)")
plt.ylabel("population")
plt.title("Species C")
plt.legend()