基于R语言非线性优化求解Walk-in中心最优人员配置问题

Got it, let's tackle this Walk-in center staffing problem step by step. We'll use R's nonlinear optimization tools paired with queueing theory to meet your constraints while minimizing total staff.

1. Problem Breakdown #

First, let's formalize what we're working with:

• Knowns: 15-hour daily demand (7:00-21:00, hourly breakdown), hourly service capacity per employee
• Constraints:
  • Customer wait time ≤ your set threshold (minutes)
  • All daily arrivals must be served before closing (21:00)
• Goal: Find the minimum number of staff (hourly) that satisfies all constraints

2. Modeling Approach #

We'll use the M/M/c queueing model here—it's standard for Walk-in services where arrivals follow a Poisson distribution and service times are exponential. The wait time formula from this model is nonlinear, which is why we need a nonlinear optimizer.

Key Definitions: #

• x_t: Number of staff scheduled for hour t (integer value, since you can't schedule half an employee)
• λ_t: Hourly arrival rate (your given demand for hour t)
• μ: Hourly service capacity per employee (e.g., 5 customers per hour per staff member)
• T: Wait time threshold (converted to hours for calculations, then back to minutes)

Objective & Constraints: #

• Objective Function: Minimize total daily staff → min(sum(x_t))
• Constraints:
  1. For every hour t: Wait time W_t ≤ T (using M/M/c wait time formula)
  2. For every hour t: x_t * μ > λ_t (prevents infinite wait times from an understaffed system)
  3. Total daily service capacity ≥ total daily demand → sum(x_t * μ) ≥ sum(λ_t) (ensures all customers are served by closing)

3. R Implementation #

We'll use the nloptr package—it supports nonlinear constraints and integer variables, which is perfect here.

Step 1: Install & Load Required Package #

install.packages("nloptr")
library(nloptr)

Step 2: Define Your Parameters #

Replace these with your actual data:

# Hourly demand (7:00-21:00, 15 values)
demand <- c(8, 12, 15, 18, 20, 17, 14, 16, 19, 22, 18, 15, 12, 9, 6)
# Hourly service capacity per employee (e.g., 5 customers/hour)
mu <- 5
# Wait time threshold (minutes)
wait_threshold <- 10

Step 3: Wait Time Calculation Function #

Implements the M/M/c total wait time (queue + service) converted to minutes:

calculate_wait_time <- function(lambda, mu, c) {
  # Handle unstable system (staff can't keep up)
  if (c * mu <= lambda) {
    return(Inf)
  }
  rho <- lambda / (c * mu)
  # Probability of 0 customers in the system
  p0 <- 1 / (sum(sapply(0:(c-1), function(k) (lambda/mu)^k / factorial(k))) + 
                (lambda/mu)^c / (factorial(c) * (1 - rho)))
  # Queue wait time (hours)
  wq <- (p0 * (lambda/mu)^c * rho) / (factorial(c) * (c * mu - lambda)^2)
  # Total wait time (queue + service) converted to minutes
  (wq + 1/mu) * 60
}

Step 4: Define Objective & Constraint Functions #

# Objective: Minimize total staff
objective <- function(x) {
  sum(x)
}

# Inequality constraints (all must be ≥ 0)
constraints <- function(x) {
  n_hours <- length(x)
  lambda <- demand
  mu_local <- mu
  threshold <- wait_threshold
  
  # Constraint 1: Wait time ≤ threshold
  wait_constraints <- sapply(1:n_hours, function(t) {
    threshold - calculate_wait_time(lambda[t], mu_local, x[t])
  })
  
  # Constraint 2: Staff can handle hourly demand (avoid unstable system)
  capacity_constraints <- sapply(1:n_hours, function(t) {
    x[t] * mu_local - lambda[t]
  })
  
  # Constraint 3: Total daily capacity ≥ total demand
  total_capacity_constraint <- sum(x)*mu_local - sum(lambda)
  
  # Combine all constraints
  c(wait_constraints, capacity_constraints, total_capacity_constraint)
}

Step 5: Run the Optimization #

We'll use a global optimization algorithm suited for integer nonlinear problems:

# Initial guess: Staff to just cover hourly demand (rounded up)
initial_x <- ceiling(demand / mu)
# Variable bounds: At least 1 staff per hour, upper limit set to double the initial guess
lower_bounds <- rep(1, length(demand))
upper_bounds <- rep(max(initial_x)*2, length(demand))

# Optimization settings
opts <- list(
  "algorithm" = "NLOPT_GN_ISRES",
  "xtol_rel" = 1e-4,
  "maxeval" = 10000,
  "print_level" = 2  # Set to 0 to suppress output
)

# Run optimization (all variables are integers)
result <- nloptr(
  x0 = initial_x,
  eval_f = objective,
  eval_g_ineq = constraints,
  lb = lower_bounds,
  ub = upper_bounds,
  opts = opts,
  integer = 1:length(initial_x)
)

Step 6: Analyze Results #

# Print optimal staffing (7:00-21:00)
cat("Optimal Hourly Staffing (7:00 to 21:00):\n")
print(result$solution)
cat("\nTotal Daily Staff Required:", sum(result$solution), "\n")

# Verify wait times for each hour
wait_times <- sapply(1:length(demand), function(t) {
  calculate_wait_time(demand[t], mu, result$solution[t])
})
cat("\nHourly Wait Times (minutes):\n")
print(round(wait_times, 2))

4. Customization Tips #

• Adjust Queue Model: If your service times don't follow an exponential distribution, swap the calculate_wait_time function for an M/G/c model (you'll need to know the variance of service times).
• Add Staff Continuity: If staff must work consecutive shifts (e.g., 4-hour blocks), add constraints linking x_t to x_{t+1}, x_{t+2}, etc.
• Speed Up Calculations: If you're short on time, run a continuous optimization first, round to integers, then verify if constraints are met (adjust manually if needed).

内容的提问来源于stack exchange，提问作者ashwini mahapatra