import matplotlib.pyplot as plt
import numpy as np

F = 0.2
V = 1
t_transp = 5.0  # Transport delay

cA_min = 0
cA_max = np.inf

dt = 0.1
t_start = 0
t_stop = 50
N_sim = int((t_stop - t_start)/dt) + 1

t_array = np.zeros(N_sim)
cA_array = np.zeros(N_sim)
wA_non_delayed_array = np.zeros(N_sim)
wA_delayed_array = np.zeros(N_sim)

# Array for transport delay:
wA_delayed_init = 0
Nf = int(round(t_transp/dt)) + 1
delay_array = np.zeros(Nf) + wA_delayed_init

cA_k = cA_init = 0

for k in range(0, N_sim):

    t_k = k*dt
    
    if (0 <= t_k <= 10): wA_non_delayed_k = 0
    else: wA_non_delayed_k = 0.8
    
    # Moving array elements one step:
    wA_delayed_k = delay_array[-1]
    delay_array[1:] = delay_array[0:-1]
    delay_array[0] = wA_non_delayed_k

    
    dcA_dt_k = (1/V)*(wA_delayed_k - cA_k*F)
    cA_kp1 = cA_k + dt*dcA_dt_k
    cA_kp1 = np.clip(cA_kp1, cA_min, cA_max)
    
    t_array[k] = t_k
    cA_array[k] = cA_k
    wA_non_delayed_array[k] = wA_non_delayed_k
    wA_delayed_array[k] = wA_delayed_k

    cA_k = cA_kp1

plt.close('all')
plt.figure(1, figsize=(12, 9))

plt.subplot(2, 1, 1)
plt.plot(t_array, cA_array, 'g', label='cA')
plt.legend()
plt.xlabel('t [min]')
plt.ylabel('[kg/m3]')
plt.grid()

plt.subplot(2, 1, 2)
plt.plot(t_array, wA_non_delayed_array, 'b',
         label='wA_non_delayed')
plt.plot(t_array, wA_delayed_array, 'r',
         label='wA_delayed')
plt.legend()
plt.xlabel('t [min]')
plt.ylabel('[kg/min]')
plt.grid()

# plt.savefig('plot_sim_blandetank_transport.pdf')
plt.show()