simulation.py

import heapq
from collections import deque
import math
from enum import Enum
import numpy as np
import matplotlib.pyplot as plt





class EventTypes(Enum):
    ARRIVAL = 0
    BROWS_COMPLETE = 1
    LEAVE = 2

class CustomerMode(Enum):
    DAILY = 0
    PURPOSEFUL = 1
    WONDERING = 2




base_browsing_t = {
    CustomerMode.DAILY: 3,
    CustomerMode.PURPOSEFUL: 4,
    CustomerMode.WONDERING: 5
}

per_item_browsing_t = {
    CustomerMode.DAILY: 0.5,
    CustomerMode.PURPOSEFUL: 1,
    CustomerMode.WONDERING: 2
}

avg_item_count = np.array([
    5, 15, 15
])

browsing_t_noise_shape = 5   # 1/Randomness
browsing_t_noise_scale = 0.4 # Avg t / (1/randomness)




base_scan_t = 1
per_item_scan_t = 0.04
scan_noise_shape = 5   # 1/Randomness
scan_noise_scale = 1/20 # Avg t / (1/randomness)



free_cashires = 0
events = []
queue = deque()

def get_cid_generator():
    id = 0
    while True:
        yield id
        id += 1



def create_arrivals(t_start, arrival_curve_x, arrival_curve_y, customer_count, customer_mode_p):
    global arrival_over_t, total_customers
    total_customers += customer_count

    customer_mode_p = np.array(customer_mode_p)
    customer_mode_p = customer_mode_p / customer_mode_p.sum()

    # arrivals = np.cumsum(np.random.exponential(1/arrivals_per_minute, customer_count)) + t_start

    arrival_curve_y = arrival_curve_y / arrival_curve_y.sum()
    arrivals = np.random.choice(arrival_curve_x, size= customer_count, p= arrival_curve_y) + t_start # Uses a custom curve

    modes = np.random.choice(range(3), size= customer_count, p= customer_mode_p)
    item_counts = np.random.poisson(avg_item_count[modes], customer_count)


    for t, mode, item_count in zip(arrivals, modes, item_counts):
        cid = next(get_cid)
        arrival_over_t[cid] = t
        heapq.heappush(events, (t, EventTypes.ARRIVAL, cid, CustomerMode(mode), item_count))

def calculate_browsing_t(event):
    t_start, _, cid, mode, item_count = event

    noise = np.random.gamma(browsing_t_noise_shape, browsing_t_noise_scale)
    t = t_start + base_browsing_t[mode] + ( per_item_browsing_t[mode] *  item_count) + noise
    heapq.heappush(events, (t, EventTypes.BROWS_COMPLETE, cid, item_count))


def calculate_scanning_t(t_start, cid, item_count):
    global free_cashires, queue_over_t, t_queue_total, queue_over_t
    t_queue_total += t_start - queue_over_t[cid][0]
    queue_over_t[cid] = (t_start, len(queue))
    free_cashires -= 1
    noise = np.random.gamma(scan_noise_shape, scan_noise_scale)
    t = t_start + base_scan_t + ( per_item_scan_t * item_count ) + noise
    heapq.heappush(events, (t, EventTypes.LEAVE, cid))


def free_cashire():
    global free_cashires
    free_cashires += 1



def do_event():
    global t, t_per_customer_total, t_queue_total, queue_len_total, queue_len_max, queue_over_t, arrival_over_t, free_cashires

    event = heapq.heappop(events)
    t, event_type, cid = event[:3]

    if event_type == EventTypes.ARRIVAL:
        calculate_browsing_t(event)

    elif event_type == EventTypes.BROWS_COMPLETE:
        queue.append((cid, event[3])) # (cid, item_count)

        queue_len = len(queue)
        if queue_len > queue_len_max:
            queue_len_max = queue_len
        queue_len_total += queue_len

        queue_over_t[cid] = (t, queue_len)
        
        if free_cashires > 0:
            
            next_customer = queue.popleft()
            calculate_scanning_t(t, next_customer[0], next_customer[1]) # (... , cid, item_count)

    elif event_type == EventTypes.LEAVE:
        free_cashire()
        if len(queue) > 0:
            next_customer = queue.popleft()
            calculate_scanning_t(t, next_customer[0], next_customer[1]) # (... , cid, item_count)
        t_per_customer_total += t - arrival_over_t[cid]

get_cid = None
t_per_customer_total = 0
t_queue_total = 0

queue_len_total = 0
queue_len_max = 0

total_queue_over_t = []
total_arrival_over_t = []

total_customers = 0

arrival_over_t = {}
queue_over_t= {}

def simulate(gen_count, cashires, resolution_arrival, gens_to_show, arrival_curve_x, arrival_curve_y, customer_count, t_start):
    global get_cid, t_per_customer_total, t_queue_total, queue_len_total, queue_len_max, total_queue_over_t, total_arrival_over_t, total_customers, arrival_over_t, queue_over_t, free_cashires
    t_per_customer_total = 0
    t_queue_total = 0

    queue_len_total = 0
    queue_len_max = 0

    total_queue_over_t = []
    total_arrival_over_t = []

    total_customers = 0
    free_cashires = cashires
    for i in range(gen_count):
        t = 0
        get_cid = get_cid_generator()
        arrival_over_t = {}
        queue_over_t= {}

        create_arrivals(t_start, arrival_curve_x, arrival_curve_y, customer_count, [1, 3, 2])
        while len(events) > 0:
            do_event()
        total_queue_over_t.append(queue_over_t)
        total_arrival_over_t.append(arrival_over_t)


    t_max = max(
        max(value[0] for value in queue_over_t.values())
        for queue_over_t in total_queue_over_t
    )

    time_grid_queue_len = np.linspace(0, t_max, int((t_max)+1))
    time_grid_arrival = np.linspace(0, t_max, int(t_max/resolution_arrival)+1)

    total_remmaped_queue_len = []
    total_remmaped_arrival = []

    fig, ax_queue_len = plt.subplots(figsize= (10, 4))
    ax_arrival = ax_queue_len.twinx()
    for i, queue_over_t, arrival_over_t in zip(range(len(total_queue_over_t)), total_queue_over_t, total_arrival_over_t):
        values = queue_over_t.values()
        times = np.array([value[0] for value in values])
        queue_lens = np.array([value[1] for value in values])
        if gens_to_show != 0 and i % (gen_count//gens_to_show) == 0:
            ax_queue_len.step(times, queue_lens, where= "post", alpha=2/(gens_to_show+4))
        idx = np.searchsorted(times, time_grid_queue_len, side= "right") - 1
        idx[idx < 0] = 0
        remapped_queue_len = queue_lens[idx]
        total_remmaped_queue_len.append(remapped_queue_len)


        times = np.array(list(arrival_over_t.values()))
        count, _ = np.histogram(times, bins=time_grid_arrival)
        if gens_to_show != 0 and i % (gen_count//gens_to_show) == 0:
            ax_arrival.step(time_grid_arrival[:-1], count, where="post", alpha=2/(gens_to_show+4))
        total_remmaped_arrival.append(count)

    total_remmaped_queue_len = np.array(total_remmaped_queue_len)
    avg_queue_over_t = total_remmaped_queue_len.mean(axis=0)
    ax_queue_len.step(time_grid_queue_len, avg_queue_over_t, where="post", linewidth=2, color= "#ff0000", label="Average Queue Len")

    total_remmaped_arrival = np.array(total_remmaped_arrival)
    avg_arrival_over_t = total_remmaped_arrival.mean(axis=0)
    ax_arrival.step(time_grid_arrival[:-1], avg_arrival_over_t, where="post", linewidth=2, color= "#f700ff", label="Average Arrivals")


    ax_queue_len.set_xlabel("Time")
    ax_queue_len.set_ylabel("Queue Len.")
    ax_arrival.set_ylabel(f"Arrivals per {resolution_arrival}m")
    ax_queue_len.set_title(f"Simulation: (cashires = {free_cashires})")

    lines_queue, labels_queue = ax_queue_len.get_legend_handles_labels()
    lines_arrival, labels_arrival = ax_arrival.get_legend_handles_labels()
    ax_queue_len.legend(lines_queue + lines_arrival, labels_queue + labels_arrival, loc="upper right")

    avg_t = t_per_customer_total/total_customers
    avg_queue_t = t_queue_total/total_customers
    avg_queue_len = int(queue_len_total/total_customers//free_cashires)
    max_queue_len = queue_len_max//free_cashires

    return fig, (avg_t, avg_queue_t, avg_queue_len, max_queue_len)