How to use SIMMER for basic queueing

[swaheera]

I was reading about the Simmer package in R. It looks very interesting!

I am working on the following simulation problem:

• Suppose there is a factory with 5 employees
• Each day, boxes arrive at the factory according to some random date
• Employees work on the boxes based on a First In First Out (FIFO) rule : Boxes that have been waiting for longer periods of time are processed before boxes waiting for shorter periods of time
• The employees can finish a box based on some other random rate. Once they finish a box, they immediately move to the next box (employees are never idle and never take breaks)
• In the situation that there more than two employees are sitting idle since there are no new boxes, the employee that was sitting idle the longest will be assigned in the newest box
• Imagine this factory is open for 10 hours each day. I want to simulate the detailed evolution of this factory over 10 days

Is it possible to do this problem in Simmer?

Thanks,
S

Here is the code I am currently using in base R:

library(data.table)

# Step 1: Define initial parameters
lambda_arrival <- 5  # Arrival rate (Poisson)
lambda_service <- 3  # Service rate (Exponential)
hours_per_day <- 10
days <- 10
workers <- 5

results <- vector("list", days * hours_per_day)
completed_boxes <- vector("list", days * hours_per_day)

box_id_counter <- 1

# Step 2: Define a Queue Function
simulate_hour <- function(hour, day, queue, workers_status) {
    # Generate arrivals
    arrivals <- rpois(1, lambda_arrival)
    arrival_times <- rep(hour, arrivals)
    arrival_days <- rep(day, arrivals)
    arrival_ids <- seq(box_id_counter, box_id_counter + arrivals - 1)
    box_id_counter <<- box_id_counter + arrivals
   
    # Update queue with new arrivals
    if (arrivals > 0) {
        new_boxes <- data.table(box_id = arrival_ids, arrival_hour = arrival_times, arrival_day = arrival_days, arrival_id = arrival_ids, waiting_time = 0)
        queue <- rbind(queue, new_boxes)
    }
   
    # Process boxes
    completed <- data.table(worker = character(0), box_id = numeric(0), completion_hour = numeric(0), completion_day = numeric(0), arrival_hour = numeric(0), arrival_day = numeric(0))
    time_remaining <- 1  # Start with 1 hour available
   
    while (time_remaining > 0 && nrow(queue) > 0) {
        for (i in 1:workers) {
            if (nrow(queue) > 0) {
                service_time <- rexp(1, lambda_service)
                if (service_time <= time_remaining) {
                    time_remaining <- time_remaining - service_time
                    completed <- rbind(completed, data.table(worker = paste0("worker_", i), box_id = queue[1, box_id], completion_hour = hour, completion_day = day, arrival_hour = queue[1, arrival_hour], arrival_day = queue[1, arrival_day]))
                    queue <- queue[-1]  # Remove the processed box from the queue
                } else {
                    # Update waiting time if the service time extends beyond the current hour
                    queue[1, waiting_time := waiting_time + time_remaining]
                    time_remaining <- 0  # No more time left in the hour
                    break
                }
            }
        }
    }
   
    # Update waiting times for remaining boxes
    queue[, waiting_time := waiting_time + 1]
   
    # Store results
    results[[hour + (day - 1) * hours_per_day]] <<- data.table(current_hour = hour, current_day = day, queue)
    completed_boxes[[hour + (day - 1) * hours_per_day]] <<- completed
   
    return(queue)
}

# Step 3: Simulate 10 days
queue <- data.table(box_id = numeric(0), arrival_hour = numeric(0), arrival_day = numeric(0), arrival_id = numeric(0), waiting_time = numeric(0))
for (day in 1:days) {
    for (hour in 1:hours_per_day) {
        queue <- simulate_hour(hour, day, queue, workers)
    }
}

# Step 4: Summarize the factory status
summarize_factory_status <- function(results, completed_boxes, hours_per_day, days) {
    summary <- data.table(
        day = rep(1:days, each = hours_per_day),
        hour = rep(1:hours_per_day, days),
        opening_inventory = 0,
        arrivals = 0,
        completed_boxes = 0,
        closing_inventory = 0
    )
   
    for (i in 1:(days * hours_per_day)) {
        current_result <- results[[i]]
        current_completed <- completed_boxes[[i]]
       
        if (i == 1) {
            summary[i, opening_inventory := 0]
        } else {
            summary[i, opening_inventory := summary[i-1, closing_inventory]]
        }
       
        summary[i, arrivals := nrow(current_result[arrival_hour == current_result[1, current_hour] &
                                                       arrival_day == current_result[1, current_day]])]
        summary[i, completed_boxes := nrow(current_completed)]
        summary[i, closing_inventory := nrow(current_result)]
    }
   
    return(summary)
}

# Generate the summary
factory_summary <- summarize_factory_status(results, completed_boxes, hours_per_day, days)

# Print the summary data frame
print(data.frame(factory_summary))

[Enchufa2]

What you are describing is an M/M/5 system, and this can be easily implemented as shown in this vignette.

[swaheera]

Hi Dr. Ucar, Thank you so much for your reply. I should have clarified in my question: I am interesting in simulating the transient behavior of the queue. I am not interested in simulating the long term behavior of the queue. This should still be possible in the SIMMER package, correct? Thanks, S

…

________________________________ From: Iñaki Ucar ***@***.***> Sent: August 21, 2024 8:05 AM To: r-simmer/simmer ***@***.***> Cc: swaheera ***@***.***>; Author ***@***.***> Subject: Re: [r-simmer/simmer] How to use SIMMER for basic queueing (Issue #309) What you are describing is an M/M/5 system, and this can be easily implemented as shown in this vignette<https://r-simmer.org/articles/simmer-06-queueing#mmck-systems>. — Reply to this email directly, view it on GitHub<#309 (comment)>, or unsubscribe<https://github.com/notifications/unsubscribe-auth/ASHV6NCICAVKGAT5EMAZLVLZSR7BXAVCNFSM6AAAAABM3IU3POVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDGMBRHA4DQNZRGE>. You are receiving this because you authored the thread.Message ID: ***@***.***>

[Enchufa2]

Not sure if I understand what you want; simmer simulates individual arrivals, not long-term behavior. Long-term behavior can be inferred from the post-processing of the simulation data. You have every arrival, every enqueueing event, every serving event, and so on, available when the simulation ends. It would help if you could skim through the various vignettes available, which will for sure clarify how simmer works, and then we can address here specific questions.

So is it possible to simulate _____ with simmer? The answer is 99.9% yes.

[Enchufa2]

Not sure if I understand what you want

First of all, because I don't understand your simulation above. I get:

> print(data.frame(factory_summary))
    day hour opening_inventory arrivals completed_boxes closing_inventory
1     1    1                 0        3               4                 3
2     1    2                 3        3               4                 3

How is that there are more completed boxes than arrivals? How is that the closing inventory doesn't change in the second line? Etc.