R의 몬테카를로 시뮬레이션 지수 랜덤 변수입니다 ( j

다음 운동을 해결하려고하지만 실제로이 작업을 시작하는 방법에 대한 실마리는 없습니다. 내 책에서 코드처럼 보이는 코드를 찾았지만 완전히 다른 연습이며 서로 연관시키는 방법을 모르겠습니다. 도착 시뮬레이션을 시작하려면 어떻게해야하며 도착이 언제인지 어떻게 알 수 있습니까? 나는 그것들을 저장하는 방법을 알고 그것에 따라 a, b, c, d를 계산합니다. 그러나 실제로 몬테 카를로 시뮬레이션을 어떻게 시뮬레이트해야하는지 모르겠습니다. 누군가 나를 도와 줄 수 있습니까? 나는 이것이 귀하의 질문에 대한 답변이 아니라 대신 해결되는 곳이라는 것을 알고 있습니다. 그러나 문제는 시작하는 방법을 모른다는 것입니다.

IT 지원 헬프 데스크는 고객으로부터 전화를받는 5 명의 비서가있는 대기열 시스템을 나타냅니다. 통화는 Poisson 프로세스에 따라 45 초마다 평균 한 번의 통화 비율로 발생합니다. 첫 번째, 두 번째, 세 번째, 네 번째 및 다섯 번째 보조자의 서비스 시간은 모두 각각 λ1 = 0.1, λ2 = 0.2, λ3 = 0.3, λ4 = 0.4 및 λ5 = 0.5 min-1의 매개 변수를 갖는 지수 랜덤 변수입니다 ( j 번째 헬프 데스크 어시스턴트는 λk = k / 10 min-1입니다). 지원을받는 고객 외에도 최대 10 명의 다른 고객을 보류 할 수 있습니다. 이 용량에 도달하면 새 발신자는 통화 중 신호를받습니다. Monte Carlo 방법을 사용하여 다음 성능 특성을 추정하십시오.

(a) 통화 중 신호를받는 고객의 비율;

(b) 예상 응답 시간;

(c) 평균 대기 시간;

(d) 각 헬프 데스크 어시스턴트가 서비스하는 고객의 일부;

편집 : 내가 지금까지 가지고있는 것은 (별로) :

pa = 1/45sec-1

jobs = rep(1,5); onHold = rep(1,10);

jobsIndex = 0;

onHoldIndex = 0;

u = runif(1)
for (i in 1:1000) {

    if(u  <= pa){ # new arrival

        if(jobsIndex < 5) # assistant is free, #give job to assistant

            jobsIndex++;

        else #add to onHold array

            onHoldIndex++;
    }
}


답변

이것은 수행 할 수 있는 가장 유익하고 재미있는 시뮬레이션 중 하나입니다. 컴퓨터에서 독립 에이전트 를 생성 하고 상호 작용하게하고 그들이하는 일을 추적하고 발생하는 일을 연구합니다. 복잡한 시스템, 특히 순전히 수학적 분석으로는 이해할 수없는 시스템에 대해 배울 수있는 놀라운 방법입니다.

이러한 시뮬레이션을 구성하는 가장 좋은 방법은 하향식 디자인입니다.

가장 높은 수준 에서 코드는 다음과 같아야합니다.

initialize(...)
while (process(get.next.event())) {}

(이 예제와 이후의 모든 예제는 의사 코드가 아닌 실행 가능한 R 코드입니다.) 루프는 이벤트 중심 시뮬레이션입니다. get.next.event()관심있는 “이벤트”를 찾아서 이에 대한 설명을 전달 process합니다. 그것에 관한 정보). TRUE일이 잘 진행되는 한 돌아옵니다 . 오류 또는 시뮬레이션의 끝을 식별 FALSE하면 루프를 종료하고를 반환 합니다.

뉴욕시에서 결혼 허가를 기다리는 사람들이나 운전 면허증 또는 기차표를 거의 모든 곳에서 기다리는 것과 같이이 대기열을 실제로 구현한다고 상상하면 고객과 “보조원”(또는 서버)의 두 가지 에이전트를 생각합니다. . 고객은 표시하여 자신을 발표합니다. 보조원은 조명을 켜거나 사인을 보거나 창을 열어 가용성을 발표합니다. 이들은 처리해야 할 두 가지 종류의 이벤트입니다.

이러한 시뮬레이션을위한 이상적인 환경은 객체를 변경할 수 있는 진정한 객체 지향 환경입니다. 객체가 주변의 사물에 독립적으로 반응하도록 상태를 변경할 수 있습니다. R이것에 대해 절대적으로 끔찍합니다 (Fortran이 더 나을 것입니다!). 그러나 약간의주의를 기울이면 여전히 사용할 수 있습니다. 비결은 모든 정보를 공통된 일련의 데이터 구조로 유지하는 것인데,이 데이터 구조는 여러 개의 상호 작용하는 개별 절차에 의해 액세스 (및 수정) 될 수 있습니다. 그런 데이터에 변수 이름 IN ALL CAPS를 사용하는 규칙을 채택하겠습니다.

하향식 디자인 의 다음 단계 는 코딩하는 것 process입니다. 단일 이벤트 설명자에 응답합니다 e.

process <- function(e) {
  if (is.null(e)) return(FALSE)
  if (e$type == "Customer") {
    i <- find.assistant(e$time)
    if (is.null(i)) put.on.hold(e$x, e$time) else serve(i, e$x, e$time)
  } else {
    release.hold(e$time)
  }
  return(TRUE)
}

get.next.event보고 할 이벤트가없는 경우 널 이벤트에 응답해야합니다 . 그렇지 않으면 process시스템의 “비즈니스 규칙”을 구현합니다. 실제로 질문의 설명에서 스스로 작성합니다. 작동 방식은 결국 서브 루틴을 코딩 put.on.hold하고 release.hold(고객 대기 대기열 serve구현 ) (고객 지원 상호 작용을 구현 해야 함)을 제외하고는 약간의 주석이 필요합니다 .

“이벤트”란 무엇입니까? 여기에는 누가 행동 하는지 , 어떤 종류의 행동을 하는지 , 언제 발생 하는지 에 대한 정보가 포함되어야합니다 . 따라서 내 코드는이 세 종류의 정보가 포함 된 목록을 사용합니다. 그러나 get.next.event시간 만 검사하면됩니다. 이벤트 큐를 유지 관리하는 것만 담당합니다.

  1. 이벤트가 수신되면 큐에 넣을 수 있습니다.

  2. 대기열에서 가장 빠른 이벤트를 쉽게 추출하여 발신자에게 전달할 수 있습니다.

이것의 가장 좋은 구현 우선 순위 큐는힙, 하지만 너무 까다로운에서입니다 R. Norman Matloff의 The Art of R Programming (보다 유연하고 추상적이지만 제한된 큐 시뮬레이터를 제공하는) 의 제안에 따라 데이터 프레임을 사용하여 이벤트를 보유하고 레코드 중에서 최소 시간 동안 간단히 검색했습니다.

get.next.event <- function() {
  if (length(EVENTS$time) <= 0) new.customer()               # Wait for a customer$
  if (length(EVENTS$time) <= 0) return(NULL)                 # Nothing's going on!$
  if (min(EVENTS$time) > next.customer.time()) new.customer()# See text
  i <- which.min(EVENTS$time)
  e <- EVENTS[i, ]; EVENTS <<- EVENTS[-i, ]
  return (e)
}

이것이 코딩 될 수있는 방법은 여러 가지가 있습니다. 여기에 표시된 최종 버전 process은 “Assistant”이벤트에 반응하는 new.customer방식 과 작동 방식을 코딩 할 때 선택한 사항을 반영 합니다. get.next.event고객을 보류 대기열에서 꺼내고 다른 이벤트를 기다립니다. 때로는 두 가지 방법으로 새로운 고객을 찾아야 할 수도 있습니다. 첫째, 고객이 문을 기다리고 있는지 (있는 그대로) 확인하고, 둘째, 우리가보고 있지 않을 때 고객이 들어 왔는지 여부를 확인하는 것입니다.

분명 new.customer하고 next.customer.time중요한 루틴 이므로 다음에 처리하도록하겠습니다.

new.customer <- function() {
  if (CUSTOMER.COUNT < dim(CUSTOMERS)[2]) {
    CUSTOMER.COUNT <<- CUSTOMER.COUNT + 1
    insert.event(CUSTOMER.COUNT, "Customer",
                 CUSTOMERS["Arrived", CUSTOMER.COUNT])
  }
  return(CUSTOMER.COUNT)
}
next.customer.time <- function() {
  if (CUSTOMER.COUNT < dim(CUSTOMERS)[2]) {
    x <- CUSTOMERS["Arrived", CUSTOMER.COUNT]
  } else {x <- Inf}
  return(x) # Time when the next customer will arrive
}

CUSTOMERS열에 각 고객에 대한 데이터가있는 2D 배열입니다. 시뮬레이션 중에 고객을 설명 하고 경험을 기록 하는 네 개의 행 (필드로 작동)이 있습니다 . “도착 됨”, “서빙 됨”, “지속 기간”및 “지원자” 그들, 그리고 다른 방법 -1으로 바쁜 신호). 매우 유연한 시뮬레이션에서 이러한 열은 동적으로 생성되지만 R작동 방식 이 좋아서 도착 시간이 이미 무작위로 생성 된 단일 대형 매트릭스로 모든 고객을 처음에 생성하는 것이 편리합니다. next.customer.time이 행렬의 다음 열을 들여다 보면 누가 다음에 올지 알 수 있습니다. 글로벌 변수CUSTOMER.COUNT도착한 마지막 고객을 나타냅니다. 고객은이 포인터를 사용하여 매우 간단하게 관리하여 새로운 고객을 확보하기 위해 발전하고 다음 고객을 들여다보기 위해이를 넘어서 (진행하지 않고) 살펴 봅니다.

serve 시뮬레이션에서 비즈니스 규칙을 구현합니다.

serve <- function(i, x, time.now) {
  #
  # Serve customer `x` with assistant `i`.
  #
  a <- ASSISTANTS[i, ]
  r <- rexp(1, a$rate)                       # Simulate the duration of service
  r <- round(r, 2)                           # (Make simple numbers)
  ASSISTANTS[i, ]$available <<- time.now + r # Update availability
  #
  # Log this successful service event for later analysis.
  #
  CUSTOMERS["Assistant", x] <<- i
  CUSTOMERS["Served", x] <<- time.now
  CUSTOMERS["Duration", x] <<- r
  #
  # Queue the moment the assistant becomes free, so they can check for
  # any customers on hold.
  #
  insert.event(i, "Assistant", time.now + r)
  if (VERBOSE) cat(time.now, ": Assistant", i, "is now serving customer",
                   x, "until", time.now + r, "\n")
  return (TRUE)
}

이것은 간단합니다. ASSISTANTS는 두 개의 필드가있는 데이터 프레임입니다 : capabilities(서비스 요금 제공) 및 available. 다음에 보조자가 해제 될 시간을 표시합니다. 어시스턴트의 기능에 따라 임의의 서비스 기간을 생성하고, 어시스턴트가 다음에 사용 가능하게되는 시간을 업데이트하고, CUSTOMERS데이터 구조 에 서비스 간격을 로깅하여 고객에게 서비스를 제공합니다 . 이 VERBOSE플래그는 테스트 및 디버깅에 편리합니다. true 인 경우 주요 처리 지점을 설명하는 일련의 영어 문장을 내 보냅니다.

비서가 고객에게 할당되는 방법은 중요하고 흥미 롭습니다. 임의의 순서, 고정 된 순서, 또는 가장 긴 (또는 가장 짧은) 시간에 누가 자유 로웠는지에 따라 여러 가지 절차를 상상할 수 있습니다. 이들 중 많은 부분이 주석 처리 된 코드로 설명되어 있습니다.

find.assistant <- function(time.now) {
  j <- which(ASSISTANTS$available <= time.now)
  #if (length(j) > 0) {
  #  i <- j[ceiling(runif(1) * length(j))]
  #} else i <- NULL                                    # Random selection
  #if (length(j) > 0) i <- j[1] else i <- NULL         # Pick first assistant
  #if (length(j) > 0) i <- j[length(j)] else i <- NULL # Pick last assistant
  if (length(j) > 0) {
    i <- j[which.min(ASSISTANTS[j, ]$available)]
  } else i <- NULL                                     # Pick most-rested assistant
  return (i)
}

나머지 시뮬레이션은 실제로R 표준 데이터 구조 (주로 보류 큐의 순환 버퍼)를 구현 하도록 설득 하는 일상적인 연습 입니다 . 글로벌로 amok를 실행하고 싶지 않기 때문에 모든 것을 단일 절차에 넣었습니다 sim. 그 논증은 문제를 설명합니다 : 시뮬레이션 할 고객의 수 ( n.events), 고객 도착율, 보조원의 능력 및 보류 대기열의 크기 (대기열을 모두 제거하기 위해 0으로 설정할 수 있음).

r <- sim(n.events=250, arrival.rate=60/45, capabilities=1:5/10, hold.queue.size=10)

CUSTOMERSR

50

250

각 고객의 경험은 도착시 원형 기호, 대기중인 대기중인 검은 색 선, 보조자와의 상호 작용 시간 (색상 및 선 종류)에 대한 컬러 선이있는 수평 타임 라인으로 표시됩니다. 조수를 구별하십시오). 이 고객 그림 아래에는 보조 직원의 경험을 보여주는 고객 그림이 있으며, 고객과 관계를 맺고 있지 않은 시간을 표시합니다. 각 활동 간격의 끝점은 세로 막대로 구분됩니다.

로 실행 verbose=TRUE하면 시뮬레이션의 텍스트 출력은 다음과 같습니다.

...
160.71 : Customer 211 put on hold at position 1
161.88 : Customer 212 put on hold at position 2
161.91 : Assistant 3 is now serving customer 213 until 163.24
161.91 : Customer 211 put on hold at position 2
162.68 : Assistant 4 is now serving customer 212 until 164.79
162.71 : Assistant 5 is now serving customer 211 until 162.9
163.51 : Assistant 5 is now serving customer 214 until 164.05
...

160

165

고객이 통화 중 신호를 받고 있음을 나타 내기 위해 특수 (빨간색) 기호를 사용하여 고객 식별자로 보류 기간을 표시하여 고객의 보류 대기 경험을 연구 할 수 있습니다.

(이 모든 플롯이이 서비스 대기열을 관리하는 사람에게 훌륭한 실시간 대시 보드를 만드는 것은 아닙니다!)

전달 된 매개 변수를 변경했을 때 얻은 플롯과 통계를 비교하는 것은 매우 흥미 롭습니다 sim. 고객이 처리하기에 너무 빨리 도착하면 어떻게됩니까? 보류 대기열을 작게 만들거나 제거하면 어떻게됩니까? 어시스턴트를 다른 방식으로 선택하면 어떤 변화가 있습니까? 보조자의 수와 기능은 고객 경험에 어떤 영향을 미칩니 까? 일부 고객이 멀어 지거나 오랫동안 보류되기 시작하는 중요한 포인트는 무엇입니까?


일반적으로 이와 같은 명백한 자기 연구 질문에 대해서는 여기서 멈추고 나머지 세부 사항을 연습으로 남겨 두십시오. 그러나, 나는 지금까지 이것을 얻었을 수도 있고 독자적으로 이것을 시도하고 (그리고 다른 목적을 위해 그것을 수정하고) 관심이있는 독자들을 실망시키지 않기를 원합니다. 그래서 아래에 전체 작업 코드가 추가됩니다.

티이자형엑스

$
sim <- function(n.events, verbose=FALSE, ...) {
  #
  # Simulate service for `n.events` customers.
  #
  # Variables global to this simulation (but local to the function):
  #
  VERBOSE <- verbose         # When TRUE, issues informative message
  ASSISTANTS <- list()       # List of assistant data structures
  CUSTOMERS <- numeric(0)    # Array of customers that arrived
  CUSTOMER.COUNT <- 0        # Number of customers processed
  EVENTS <- list()           # Dynamic event queue   
  HOLD <- list()             # Customer on-hold queue
  #............................................................................#
  #
  # Start.
  #
  initialize <- function(arrival.rate, capabilities, hold.queue.size) {
    #
    # Create common data structures.
    #
    ASSISTANTS <<- data.frame(rate=capabilities,     # Service rate
                              available=0            # Next available time
    )
    CUSTOMERS <<- matrix(NA, nrow=4, ncol=n.events,
                         dimnames=list(c("Arrived",  # Time arrived
                                         "Served",   # Time served
                                         "Duration", # Duration of service
                                         "Assistant" # Assistant id
                         )))
    EVENTS <<- data.frame(x=integer(0),              # Assistant or customer id
                          type=character(0),         # Assistant or customer
                          time=numeric(0)            # Start of event
    )
    HOLD <<- list(first=1,                           # Index of first in queue
                  last=1,                            # Next available slot
                  customers=rep(NA, hold.queue.size+1))
    #
    # Generate all customer arrival times in advance.
    #
    CUSTOMERS["Arrived", ] <<- cumsum(round(rexp(n.events, arrival.rate), 2))
    CUSTOMER.COUNT <<- 0
    if (VERBOSE) cat("Started.\n")
    return(TRUE)
  }
  #............................................................................#
  #
  # Dispatching.
  #
  # Argument `e` represents an event, consisting of an assistant/customer 
  # identifier `x`, an event type `type`, and its time of occurrence `time`.
  #
  # Depending on the event, a customer is either served or an attempt is made
  # to put them on hold.
  #
  # Returns TRUE until no more events occur.
  #
  process <- function(e) {
    if (is.null(e)) return(FALSE)
    if (e$type == "Customer") {
      i <- find.assistant(e$time)
      if (is.null(i)) put.on.hold(e$x, e$time) else serve(i, e$x, e$time)
    } else {
      release.hold(e$time)
    }
    return(TRUE)
  }#$
  #............................................................................#
  #
  # Event queuing.
  #
  get.next.event <- function() {
    if (length(EVENTS$time) <= 0) new.customer()
    if (length(EVENTS$time) <= 0) return(NULL)
    if (min(EVENTS$time) > next.customer.time()) new.customer()
    i <- which.min(EVENTS$time)
    e <- EVENTS[i, ]; EVENTS <<- EVENTS[-i, ]
    return (e)
  }
  insert.event <- function(x, type, time.occurs) {
    EVENTS <<- rbind(EVENTS, data.frame(x=x, type=type, time=time.occurs))
    return (NULL)
  }
  # 
  # Customer arrivals (called by `get.next.event`).
  #
  # Updates the customers pointer `CUSTOMER.COUNT` and returns the customer
  # it newly points to.
  #
  new.customer <- function() {
    if (CUSTOMER.COUNT < dim(CUSTOMERS)[2]) {
      CUSTOMER.COUNT <<- CUSTOMER.COUNT + 1
      insert.event(CUSTOMER.COUNT, "Customer",
                   CUSTOMERS["Arrived", CUSTOMER.COUNT])
    }
    return(CUSTOMER.COUNT)
  }
  next.customer.time <- function() {
    if (CUSTOMER.COUNT < dim(CUSTOMERS)[2]) {
      x <- CUSTOMERS["Arrived", CUSTOMER.COUNT]
    } else {x <- Inf}
    return(x) # Time when the next customer will arrive
  }
  #............................................................................#
  #
  # Service.
  #
  find.assistant <- function(time.now) {
    #
    # Select among available assistants.
    #
    j <- which(ASSISTANTS$available <= time.now)
    #if (length(j) > 0) {
    #  i <- j[ceiling(runif(1) * length(j))]
    #} else i <- NULL                                    # Random selection
    #if (length(j) > 0) i <- j[1] else i <- NULL         # Pick first assistant
    #if (length(j) > 0) i <- j[length(j)] else i <- NULL # Pick last assistant
    if (length(j) > 0) {
      i <- j[which.min(ASSISTANTS[j, ]$available)]
    } else i <- NULL # Pick most-rested assistant
    return (i)
  }#$
  serve <- function(i, x, time.now) {
    #
    # Serve customer `x` with assistant `i`.
    #
    a <- ASSISTANTS[i, ]
    r <- rexp(1, a$rate)                       # Simulate the duration of service
    r <- round(r, 2)                           # (Make simple numbers)
    ASSISTANTS[i, ]$available <<- time.now + r # Update availability
    #
    # Log this successful service event for later analysis.
    #
    CUSTOMERS["Assistant", x] <<- i
    CUSTOMERS["Served", x] <<- time.now
    CUSTOMERS["Duration", x] <<- r
    #
    # Queue the moment the assistant becomes free, so they can check for
    # any customers on hold.
    #
    insert.event(i, "Assistant", time.now + r)
    if (VERBOSE) cat(time.now, ": Assistant", i, "is now serving customer",
                     x, "until", time.now + r, "\n")
    return (TRUE)
  }
  #............................................................................#
  #
  # The on-hold queue.
  #
  # This is a cicular buffer implemented by an array and two pointers,
  # one to its head and the other to the next available slot.
  #
  put.on.hold <- function(x, time.now) {
    #
    # Try to put customer `x` on hold.
    #
    if (length(HOLD$customers) < 1 ||
          (HOLD$first - HOLD$last %% length(HOLD$customers) == 1)) {
      # Hold queue is full, alas.  Log this occurrence for later analysis.
      CUSTOMERS["Assistant", x] <<- -1 # Busy signal
      if (VERBOSE) cat(time.now, ": Customer", x, "got a busy signal.\n")
      return(FALSE)
    }
    #
    # Add the customer to the hold queue.
    #
    HOLD$customers[HOLD$last] <<- x
    HOLD$last <<- HOLD$last %% length(HOLD$customers) + 1
    if (VERBOSE) cat(time.now, ": Customer", x, "put on hold at position",
                 (HOLD$last - HOLD$first - 1) %% length(HOLD$customers) + 1, "\n")
    return (TRUE)
  }
  release.hold <- function(time.now) {
    #
    # Pick up the next customer from the hold queue and place them into
    # the event queue.
    #
    if (HOLD$first != HOLD$last) {
      x <- HOLD$customers[HOLD$first]   # Take the first customer
      HOLD$customers[HOLD$first] <<- NA # Update the hold queue
      HOLD$first <<- HOLD$first %% length(HOLD$customers) + 1
      insert.event(x, "Customer", time.now)
    }
  }$
  #............................................................................#
  #
  # Summaries.
  #
  # The CUSTOMERS array contains full information about the customer experiences:
  # when they arrived, when they were served, how long the service took, and
  # which assistant served them.
  #
  summarize <- function() return (list(c=CUSTOMERS, a=ASSISTANTS, e=EVENTS,
                                       h=HOLD))
  #............................................................................#
  #
  # The main event loop.
  #
  initialize(...)
  while (process(get.next.event())) {}
  #
  # Return the results.
  #
  return (summarize())
}
#------------------------------------------------------------------------------#
#
# Specify and run a simulation.
#
set.seed(17)
n.skip <- 200  # Number of initial events to skip in subsequent summaries
system.time({
  r <- sim(n.events=50+n.skip, verbose=TRUE,
           arrival.rate=60/45, capabilities=1:5/10, hold.queue.size=10)
})
#------------------------------------------------------------------------------#
#
# Post processing.
#
# Skip the initial phase before equilibrium.
#
results <- r$c
ids <- (n.skip+1):(dim(results)[2])
arrived <- results["Arrived", ]
served <- results["Served", ]
duration <- results["Duration", ]
assistant <- results["Assistant", ]
assistant[is.na(assistant)] <- 0   # Was on hold forever
ended <- served + duration
#
# A detailed plot of customer experiences.
#
n.events <- length(ids)
n.assistants <- max(assistant, na.rm=TRUE)
colors <- rainbow(n.assistants + 2)
assistant.color <- colors[assistant + 2]
x.max <- max(results["Served", ids] + results["Duration", ids], na.rm=TRUE)
x.min <- max(min(results["Arrived", ids], na.rm=TRUE) - 2, 0)
#
# Lay out the graphics.
#
layout(matrix(c(1,1,2,2), 2, 2, byrow=TRUE), heights=c(2,1))
#
# Set up the customers plot.
#
plot(c(x.min, x.max), range(ids), type="n",
     xlab="Time", ylab="Customer Id", main="Customers")
#
# Place points at customer arrival times.
#
points(arrived[ids], ids, pch=21, bg=assistant.color[ids], col="#00000070")
#
# Show wait times on hold.
#
invisible(sapply(ids, function(i) {
  if (!is.na(served[i])) lines(x=c(arrived[i], served[i]), y=c(i,i))
}))
#
# More clearly show customers getting a busy signal.
#
ids.not.served <- ids[is.na(served[ids])]
ids.served <- ids[!is.na(served[ids])]
points(arrived[ids.not.served], ids.not.served, pch=4, cex=1.2)
#
# Show times of service, colored by assistant id.
#
invisible(sapply(ids.served, function(i) {
  lines(x=c(served[i], ended[i]), y=c(i,i), col=assistant.color[i], lty=assistant[i])
}))
#
# Plot the histories of the assistants.
#
plot(c(x.min, x.max), c(1, n.assistants)+c(-1,1)/2, type="n", bty="n",
     xlab="", ylab="Assistant Id", main="Assistants")
abline(h=1:n.assistants, col="#808080", lwd=1)
invisible(sapply(1:(dim(results)[2]), function(i) {
  a <- assistant[i]
  if (a > 0) {
    lines(x=c(served[i], ended[i]), y=c(a, a), lwd=3, col=colors[a+2])
    points(x=c(served[i], ended[i]), y=c(a, a), pch="|", col=colors[a+2])
  }
}))
#
# Plot the customer waiting statistics.
#
par(mfrow=c(1,1))
i <- is.na(served)
plot(served - arrived, xlab="Customer Id", ylab="Minutes",
     main="Service Wait Durations")
lines(served - arrived, col="Gray")
points(which(i), rep(0, sum(i)), pch=16, col="Red")
#
# Summary statistics.
#
mean(!is.na(served)) # Proportion of customers served
table(assistant)