R/RcppExports.R

In rchan26/layeredBB: Simulates layered Brownian bridges using C++ and Rcpp

# Generated by using Rcpp::compileAttributes() -> do not edit by hand
# Generator token: 10BE3573-1514-4C36-9D1C-5A225CD40393

#' Bessel Layer simulation
#'
#' Simulates a Bessel layer l for a given sequence a
#'
#' @param x start value of Brownian bridge
#' @param y end value of Brownian bridge
#' @param s start time of Brownian bridge
#' @param t end time of Brownian bridge
#' @param a vector/sequence of numbers
#'
#' @examples
#' bessel_layer_simulation(x = 0,
#'                         y = 0,
#'                         s = 0,
#'                         t = 1,
#'                         mult = 0.5)
#' 
#' @return 
#' A list with the following items:
#' \describe{
#'   \item{L}{Hard lower bound}
#'   \item{l}{Soft lower bound}
#'   \item{u}{Soft upper bound}
#'   \item{U}{Hard upper bound}
#' }
#' where the Bessel layer is [L, U] and either the minimum occurs in [L, l] or
#' the maximum occurs in [u, U]
bessel_layer_simulation <- function(x, y, s, t, mult = 1) {
    .Call(`_layeredBB_bessel_layer_simulation`, x, y, s, t, mult)
}

#' Multiple Bessel Layer simulation
#'
#' Simulates a Bessel layer l for a given sequence a for each component of the Brownian bridge
#'
#' @param dim dimension of Brownian bridge
#' @param x vector of start values of Brownian bridge (length of vector must be equal to dim)
#' @param y vector of end values of Brownian bridge (length of vector must be equal to dim)
#' @param s start time of Brownian bridge
#' @param t end time of Brownian bridge
#' @param a vector/sequence of numbers
#'
#' @examples
#' # simulate layer information for two-dimensional Brownian bridge starting 
#' # and ending at (0,0) in time [0,1]
#' multi_bessel_layer_simulation(dim = 2,
#'                               x = c(0, 0),
#'                               y = c(0, 0),
#'                               s = 0,
#'                               t = 1,
#'                               mult = 0.5)
#' 
#' @return 
#' A list of length dim where list[i] is the Bessel layer for component i,
#' which is represented in a list with the following items:
#' \describe{
#'   \item{L}{Hard lower bound}
#'   \item{l}{Soft lower bound}
#'   \item{u}{Soft upper bound}
#'   \item{U}{Hard upper bound}
#' }
#' where the Bessel layer for component i is [L, U] and either the minimum 
#' occurs in [L, l] or the maximum occurs in [u, U] 
multi_bessel_layer_simulation <- function(dim, x, y, s, t, mult = 1) {
    .Call(`_layeredBB_multi_bessel_layer_simulation`, dim, x, y, s, t, mult)
}

#' Brownian Motion path sampler 
#'
#' Simulation of a path of a Brownian motion at given times
#'
#' @param x start value of Brownian motion
#' @param times vector of real numbers to simulate Brownian motion
#'              (must be a vector of at least length two, where times[1]
#'               is the start time of the Brownian motion, and times[2:length(times)]
#'               are the times you wish to further simulate)
#'
#' @return Matrix of the simulated Brownian motion path at all
#'         included time points. The times are sorted. 
#'         The first row are the points of the Brownian motion (named 'X')
#'         second row are corresponding times (named 'times')
#'
#' @examples
#' # simulating path for Brownian motion starting at 0 between [0,1]
#' path <- Brownian_motion_path_sampler(x = 0, times = seq(0, 1, 0.01))
#' plot(x = path['time',], y = path['X',], pch = 20, xlab = 'Time', ylab = 'X')
#' lines(x = path['time',], y = path['X',])
#' 
#' # comparing the simulated distribution of simulated points to the
#' # theoretical distribution of simulated points
#' # set variables
#' x <- 0
#' start_time <- 1.8
#' end_time <- 5
#' replicates <- 10000
#' paths <- list()
#' # repeatedly simulate Brownian bridge 
#' for (i in 1:replicates) {
#'   paths[[i]] <- Brownian_motion_path_sampler(x = x, times = seq(start_time, end_time, 0.01))
#' }
#' # select the points at the specified time q
#' index <- which(seq(start_time, end_time, 0.01)==end_time)
#' simulated_points <- sapply(1:replicates, function(i) paths[[i]]['X', index])
#' # calculate the theoretical mean and standard deviation of the simulated points at time q
#' theoretical_mean <- x
#' theoretical_sd <- sqrt(end_time-start_time)
#' # plot distribution of the simulated points and the theoretical distribution
#' plot(density(simulated_points))
#' curve(dnorm(x, theoretical_mean, theoretical_sd), add = T, col = 'red')
#' print(paste('Theoretical variance is', end_time-start_time, 'and sample variance is', var(simulated_points)))
Brownian_motion_path_sampler <- function(x, times) {
    .Call(`_layeredBB_Brownian_motion_path_sampler`, x, times)
}

#' Multi-dimensional Brownian Motion path sampler
#'
#' Simulation of a multi-dimensional Brownian Motion, at given times
#'
#' @param dim dimension of Brownian motion
#' @param x start value of Brownian motion
#' @param times vector of real numbers to simulate Brownian motion
#'              (must be a vector of at least length two, where times[1]
#'               is the start time of the Brownian motion, and times[2:length(times)]
#'               are the times you wish to further simulate)
#' 
#' @return Matrix of the simulated Brownian motion path at all
#'         included time points. The times are sorted. 
#'         The first dim rows are the points of the Brownian motion
#'         dim+1 row are corresponding times
#'
#' @examples
#' # simulate two-dimensional Brownian bridge starting and ending 
#' # at (0,0) in time [0,1]
#' multi_brownian_motion(dim = 2,
#'                       x = c(0,0),
#'                       times = c(0.1, 0.2, 0.4, 0.6, 0.8))
#'                       
#' # note that the times are sorted and duplicates are removed
#' multi_brownian_motion(dim = 2,
#'                       x = c(0.5,1),
#'                       times = c(0.1, 0.2, 0.4, 0.6, 0.6, 0.8, 0.1))
multi_brownian_motion <- function(dim, x, times) {
    .Call(`_layeredBB_multi_brownian_motion`, dim, x, times)
}

#' Brownian Bridge path sampler (Algorithm 13 in ST329)
#'
#' Simulation of a path of a Brownian bridge at given times
#'
#' @param x start value of Brownian bridge
#' @param y end value of Brownian bridge
#' @param s start time of Brownian bridge
#' @param t end time of Brownian bridge
#' @param times vector of real numbers to simulate Brownian bridge
#'
#' @return A list with the following components
#' \describe{
#'   \item{full_path}{Matrix of the simulated Brownian bridge path at all 
#'                    included time points, i.e. s, t and times. The times
#'                    are sorted and duplicates are removed. The first row
#'                    are the points of the Brownian bridge (named 'X') 
#'                    second row are corresponding times (named 'time')}
#'   \item{simulated_path}{Matrix of the simulated Brownian bridge path only at 
#'                         the specified times passed into the function, i.e. 
#'                         the times vector. The times are not sorted and
#'                         duplicates are not removed. The first row
#'                         are the points of the Brownian bridge (named 'X') 
#'                         second row are corresponding times (named 'time')}
#' }
#'
#' @examples 
#' # simulating paths for time [0,1] and plotting them
#' start <- runif(1, -1, 1)
#' end <- runif(1, -1, 1)
#' path <- Brownian_bridge_path_sampler(x = start,
#'                                      y = end,
#'                                      s = 0,
#'                                      t = 1,
#'                                      times = seq(0, 1, 0.01))$full_path
#' plot(x = path['time',], y = path['X',], pch = 20, xlab = 'Time', ylab = 'X')
#' lines(x = path['time',], y = path['X',])
#' 
#' # notice that simulated_path only includes points that are included in times vector
#' # note that simulated_path does not remove duplicates passed into times
#' Brownian_bridge_path_sampler(x = 0,
#'                              y = 1,
#'                              s = 0,
#'                              t = 1,
#'                              times = c(0.1, 0.2, 0.4, 0.6, 0.6, 0.8, 0.1))
#'
#' # comparing the simulated distribution of simulated points to the
#' # theoretical distribution of simulated points
#' # set variables
#' x <- 0.53
#' y <- 4.32
#' s <- 0.53
#' t <- 2.91
#' q <- 1.72
#' replicates <- 10000
#' paths <- list()
#' # repeatedly simulate Brownian bridge 
#' for (i in 1:replicates) {
#'   paths[[i]] <- Brownian_bridge_path_sampler(x = x,
#'                                              y = y,
#'                                              s = s,
#'                                              t = t,
#'                                              times = seq(s, t, 0.01))
#' }
#' # select the points at the specified time q
#' index <- which(seq(s, t, 0.01)==q)
#' simulated_points <- sapply(1:replicates, function(i) paths[[i]]$full_path['X', index])
#' # calculate the theoretical mean and standard deviation of the simulated points at time q
#' theoretical_mean <- x + (q-s)*(y-x)/(t-s)
#' theoretical_sd <- sqrt((t-q)*(q-s)/(t-s))
#' # plot distribution of the simulated points and the theoretical distribution
#' plot(density(simulated_points))
#' curve(dnorm(x, theoretical_mean, theoretical_sd), add = T, col = 'red')
Brownian_bridge_path_sampler <- function(x, y, s, t, times) {
    .Call(`_layeredBB_Brownian_bridge_path_sampler`, x, y, s, t, times)
}

#' Multi-dimensional Brownian Bridge path sampler
#'
#' Simulation of a multi-dimensional Brownian Bridge, at given times
#'
#' @param dim dimension of Brownian bridge
#' @param x start value of Brownian bridge
#' @param y end value of Brownian bridge
#' @param s start time of Brownian bridge
#' @param t end time of Brownian bridge
#' @param times vector of real numbers to simulate Brownian bridge
#' 
#' @return A list with the following components
#' \describe{
#'   \item{full_path}{Matrix of the simulated Brownian bridge path at all 
#'                    included time points, i.e. s, t and times. The times
#'                    are sorted and duplicates are removed. The first dim rows
#'                    are the points of the Brownian bridge in each component,
#'                    last row gives the corresponding times}
#'   \item{simulated_path}{Matrix of the simulated Brownian bridge path only at 
#'                         the specified times passed into the function, i.e. 
#'                         the times vector. The times are not sorted and
#'                         duplicates are not removed. The first dim rows
#'                         are the points of the Brownian bridge in each component,
#'                         last row gives the corresponding times}
#' }
#'
#' @examples
#' # simulate two-dimensional Brownian bridge starting and ending 
#' # at (0,0) in time [0,1]
#' multi_brownian_bridge(dim = 2,
#'                       x = c(0,0),
#'                       y = c(0,0),
#'                       s = 0,
#'                       t = 1,
#'                       times = c(0.1, 0.2, 0.4, 0.6, 0.8))
#'                       
#' # note that simulated_path does not remove duplicates passed into times
#' multi_brownian_bridge(dim = 2,
#'                       x = c(0,0),
#'                       y = c(0,0),
#'                       s = 0,
#'                       t = 1,
#'                       times = c(0.1, 0.2, 0.4, 0.6, 0.6, 0.8, 0.1))
multi_brownian_bridge <- function(dim, x, y, s, t, times) {
    .Call(`_layeredBB_multi_brownian_bridge`, dim, x, y, s, t, times)
}

#' Brownian Bridge minimum point sampler (Algorithm 14 in ST329)
#'
#' Simulation of a minimum point of a Brownian bridge
#'
#' @param x start value of Brownian bridge
#' @param y end value of Brownian bridge
#' @param s start time of Brownian bridge
#' @param t end time of Brownian bridge
#' @param low_bound Lower bound of minimum point (low_bound < up_bound <= min(x,y))
#' @param up_bound Upper bound of minimum point (low_bound < up_bound <= min(x,y))
#' @param checks logical value to determine if arguments that are passed into 
#'        the function are checked. Things that are checked include that 
#'        low_bound < up_bound <= min(x,y) and that s < t
#'
#' @return vector: the simulated minimum, 'min', and time where minimum occurs, 'tau'
#'
#' @examples
#' # simulate a minimum between 0 and 1 of a Brownian bridge starting 
#' # at 0 and ending at 0 in time [0,1]
#' min_sampler(x = 0, y = 0, s = 0, t = 1, low_bound = -1, up_bound = 0)
#'
#' # plotting multiple simulated minimums and their times
#' minimums <- sapply(1:5000, function(i) {
#'   min_sampler(x = 0, y = 0, s = 0, t = 1, low_bound = -10, up_bound = 0)
#' })
#' plot(x = minimums[2,], y = minimums[1,], pch = 20, lwd = 0.1,
#'      xlab = 'Time', ylab = 'X')
min_sampler <- function(x, y, s, t, low_bound, up_bound, checks = TRUE) {
    .Call(`_layeredBB_min_sampler`, x, y, s, t, low_bound, up_bound, checks)
}

#' Bessel Bridge point sampler given minimum (Algorithm 15 in ST329)
#'
#' Simulation of a point of a Bessel bridge at time q, given minimum occurs at time tau
#'
#' @param x start value of Bessel bridge
#' @param y end value of Bessel bridge
#' @param s start time of Bessel bridge
#' @param t end time of Bessel bridge
#' @param m minimum point
#' @param tau time of minimum point
#' @param q time of simulation
#' @param checks logical value to determine if arguments that are passed into 
#'        the function are checked. Things that are checked include that
#'        s < t, that q is in [s,t], that tau is in [s,t], that m <= min(x,y) 
#'        and that if tau == s or tau == t, then m == x or m == y, respectively
#'
#' @return simulated point of the Bessel bridge at time q
#'
#' @examples
#' # simulating a point at q=0.2 for a Bessel bridge starting at 0 and ending 
#' # at 0 in time [0,1] given minimum is at -0.4 at time 0.6
#' min_Bessel_bridge_sampler(x = 0,
#'                           y = 0,
#'                           s = 0,
#'                           t = 1,
#'                           m = -0.4,
#'                           tau = 0.6,
#'                           q = 0.2)
min_Bessel_bridge_sampler <- function(x, y, s, t, m, tau, q, checks = TRUE) {
    .Call(`_layeredBB_min_Bessel_bridge_sampler`, x, y, s, t, m, tau, q, checks)
}

#' Bessel Bridge path sampler given minimum
#'
#' Simulation of a path of a Bessel bridge at given times, given minimum occurs at time tau
#'
#' @param x start value of Bessel bridge
#' @param y end value of Bessel bridge
#' @param s start time of Bessel bridge
#' @param t end time of Bessel bridge
#' @param m minimum point
#' @param tau time of minimum point
#' @param times vector of real numbers to simulate Bessel bridge
#' @param checks logical value to determine if arguments that are passed into 
#'        the function are checked. Things that are checked include that 
#'        s < t, that requested simulation times are in [s,t], that m <= min(x,y)
#'        and that if tau == s or tau == t, then m == x or m == y, respectively
#'
#' @return A list with the following components
#' \describe{
#'   \item{full_path}{Matrix of the simulated Bessel bridge path at all 
#'                    included time points, i.e. s, t and times. The times
#'                    are sorted and duplicates are removed. The first row
#'                    are the points of the Brownian bridge (named 'X') 
#'                    second row are corresponding times (named 'time')}
#'   \item{simulated_path}{Matrix of the simulated Bessel bridge path only at 
#'                         the specified times passed into the function, i.e. 
#'                         the times vector. The times are not sorted and
#'                         duplicates are not removed. The first row
#'                         are the points of the Bessel bridge (named 'X') 
#'                         second row are corresponding times (named 'time')}
#'   \item{remove_m_path}{Matrix of the simulated Bessel bridge path only at 
#'                        all included times points excluding tau. These times
#'                        are sorted and duplicates are removed. The first row
#'                        are the points of the Bessel bridge (named 'X') 
#'                        second row are corresponding times (named 'time'). 
#'                        Note that the minimum point is included if it is 
#'                        passed into the times vector}
#' }
#'
#' @examples
#' # simulating a path at times=c(0.2, 0.4, 0.8) for a Bessel bridge starting 
#' # at 0 and ending at 0 in time [0,1] given minimum is at -0.4 at time 0.6
#' min_Bessel_bridge_path_sampler(x = 0,
#'                                y = 0,
#'                                s = 0,
#'                                t = 1,
#'                                m = -0.4,
#'                                tau = 0.6,
#'                                times = c(0.2, 0.4, 0.8))
#' 
#' # note that remove_m_path will still include the minimum if passed into times
#' # also note that simulated_path does not remove duplicates passed into times
#' min_Bessel_bridge_path_sampler(x = 0,
#'                                y = 0,
#'                                s = 0,
#'                                t = 1,
#'                                m = -0.4,
#'                                tau = 0.6,
#'                                times = c(0.2, 0.4, 0.6, 0.8, 0.6))
#' 
#' # another example
#' start <- runif(1, -1, 1)
#' end <- runif(1, -1, 1)
#' min <- min_sampler(x = start,
#'                    y = end,
#'                    s = 0,
#'                    t = 1,
#'                    low_bound = min(start, end)-0.4,
#'                    up_bound = min(start, end)-0.2)
#' path <- min_Bessel_bridge_path_sampler(x = start,
#'                                        y = end,
#'                                        s = 0,
#'                                        t = 1,
#'                                        m = min['min'],
#'                                        tau = min['tau'],
#'                                        times = seq(0, 1, 0.01))$full_path
#' plot(x = path['time',], y = path['X',], pch = 20, xlab = 'Time', ylab = 'X')
#' lines(x = path['time',], y = path['X',])
#' points(x = min['tau'], y = min['min'], col = 'red', pch = 20)
min_Bessel_bridge_path_sampler <- function(x, y, s, t, m, tau, times, checks = TRUE) {
    .Call(`_layeredBB_min_Bessel_bridge_path_sampler`, x, y, s, t, m, tau, times, checks)
}

#' Brownian Bridge maximum point sampler
#'
#' Simulation of a maximum point of a Brownian bridge
#'
#' @param x start value of Brownian bridge
#' @param y end value of Brownian bridge
#' @param s start time of Brownian bridge
#' @param t end time of Brownian bridge
#' @param low_bound Lower bound of maximum point (max(x,y) <= low_bound < up_bound)
#' @param up_bound Upper bound of maximum point (max(x,y) <= low_bound < up_bound)
#' @param checks logical value to determine if arguments that are passed into 
#'        the function are checked. Things that are checked include that 
#'        max(x,y) <= low_bound < up_bound and that s < t
#'
#' @return vector: the simulated maximum, 'max', and time where maximum occurs, 'tau'
#'
#' @examples
#' # simulate a maximum between 0 and 1 of a Brownian bridge starting at 
#' # 0 and ending at 0 in time [0,1]
#' max_sampler(x = 0, y = 0, s = 0, t = 1, low_bound = 0, up_bound = 1)
#'
#' # plotting multiple simulated maximums and their times
#' maximums <- sapply(1:5000, function(i) {
#'   max_sampler(x = 0, y = 0, s = 0, t = 1, low_bound = 0 , up_bound = 10)
#' })
#' plot(x = maximums[2,], y = maximums[1,], pch = 20, lwd = 0.1,
#'      xlab = 'Time', ylab = 'X')
max_sampler <- function(x, y, s, t, low_bound, up_bound, checks = TRUE) {
    .Call(`_layeredBB_max_sampler`, x, y, s, t, low_bound, up_bound, checks)
}

#' Bessel Bridge point sampler given maximum
#'
#' Simulation of a point of a Bessel bridge at time q, given maximum occurs at time tuu
#'
#' @param x start value of Bessel bridge
#' @param y end value of Bessel bridge
#' @param s start time of Bessel bridge
#' @param t end time of Bessel bridge
#' @param m maximum point 
#' @param tau time of maximum point
#' @param q time of simulation
#' @param checks logical value to determine if arguments that are passed into 
#'        the function are checked. Things that are checked include that
#'        s < t, that q is in [s,t], that tau is in [s,t], that m >= min(x,y) 
#'        and that if tau == s or tau == t, then m == x or m == y, respectively
#' 
#' @return simulated point of the Bessel bridge at time q
#'
#' @examples
#' # simulating a point at q=0.2 for a Bessel bridge starting at 0 and ending 
#' # at 0 in time [0,1] given maximum is at 0.4 at time 0.6
#' max_Bessel_bridge_sampler(x = 0,
#'                           y = 0,
#'                           s = 0,
#'                           t = 1,
#'                           m = 0.4,
#'                           tau = 0.6,
#'                           q = 0.2)
max_Bessel_bridge_sampler <- function(x, y, s, t, m, tau, q, checks = TRUE) {
    .Call(`_layeredBB_max_Bessel_bridge_sampler`, x, y, s, t, m, tau, q, checks)
}

#' Bessel Bridge path sampler given maximum
#'
#' Simulation of a path of a Bessel bridge at given times, given maximum occurs at time tuu
#'
#' @param x start value of Bessel bridge
#' @param y end value of Bessel bridge
#' @param s start time of Bessel bridge
#' @param t end time of Bessel bridge
#' @param m maximum point 
#' @param tau time of maximum point
#' @param times vector of real numbers to simulate Bessel bridge
#' @param checks logical value to determine if arguments that are passed into 
#'        the function are checked. Things that are checked include that 
#'        s < t, that requested simulation times are in [s,t], that m >= max(x,y)
#'        and that if tau == s or tau == t, then m == x or m == y, respectivelys
#'
#' @return A list with the following components
#' \describe{
#'   \item{full_path}{Matrix of the simulated Bessel bridge path at all 
#'                    included time points, i.e. s, t and times. The times
#'                    are sorted and duplicates are removed. The first row
#'                    are the points of the Brownian bridge (named 'X') 
#'                    second row are corresponding times (named 'time')}
#'   \item{simulated_path}{Matrix of the simulated Bessel bridge path only at 
#'                         the specified times passed into the function, i.e. 
#'                         the times vector. The times are not sorted and
#'                         duplicates are not removed. The first row
#'                         are the points of the Bessel bridge (named 'X') 
#'                         second row are corresponding times (named 'time')}
#'   \item{remove_m_path}{Matrix of the simulated Bessel bridge path only at 
#'                        all included times points excluding tau. These times
#'                        are sorted and duplicates are removed. The first row
#'                        are the points of the Bessel bridge (named 'X') 
#'                        second row are corresponding times (named 'time'). 
#'                        Note that the maximum point is included if it is 
#'                        passed into the times vector}
#' }
#'
#' @examples
#' # simulating a path at times = c(0.2, 0.4, 0.8) for a Bessel bridge starting 
#' # at 0 and ending at 0 in time [0,1] given maximum is at 0.4 at time 0.6
#' max_Bessel_bridge_path_sampler(x = 0,
#'                                y = 0,
#'                                s = 0,
#'                                t = 1,
#'                                m = 0.4,
#'                                tau = 0.6,
#'                                times = c(0.2, 0.4, 0.8))
#'
#' # note that remove_m_path will still include the minimum if passed into times
#' # also note that simulated_path does not remove duplicates passed into times
#' max_Bessel_bridge_path_sampler(x = 0,
#'                                y = 0,
#'                                s = 0,
#'                                t = 1,
#'                                m = 0.4,
#'                                tau = 0.6,
#'                                times = c(0.2, 0.4, 0.6, 0.8, 0.6))
#'
#' # another example
#' start <- runif(1, -1, 1)
#' end <- runif(1, -1, 1)
#' max <- max_sampler(x = start,
#'                    y = end,
#'                    s = 0,
#'                    t = 1,
#'                    low_bound = max(start, end)+0.2,
#'                    up_bound = max(start, end)+0.4)
#' path <- max_Bessel_bridge_path_sampler(x = start,
#'                                        y = end,
#'                                        s = 0,
#'                                        t = 1,
#'                                        m = max['max'],
#'                                        tau = max['tau'],
#'                                        times = seq(0, 1, 0.01))$full_path
#' plot(x = path['time',], y = path['X',], pch = 20, xlab = 'Time', ylab = 'X')
#' lines(x = path['time',], y = path['X',])
#' points(x = max['tau'], y = max['max'], col = 'red', pch = 20)
max_Bessel_bridge_path_sampler <- function(x, y, s, t, m, tau, times, checks = TRUE) {
    .Call(`_layeredBB_max_Bessel_bridge_path_sampler`, x, y, s, t, m, tau, times, checks)
}

#' Sigma_Bar (Equation 141 in ST329)
#'
#' This function evaluates the sigma_bar function used to simulate Bessel Layers in the infinite sums
#'
#' @param j real value
#' @param x start value of Brownian bridge
#' @param y end value of Brownian bridge
#' @param s start value of Brownian bridge
#' @param t end value of Brownian bridge
#' @param l lower bound of Brownian bridge
#' @param v upper bound of Brownian bridge
#' 
#' @return real value: easigma_bar evaluated at point j
#'
#' @examples
#' easigma_bar(j = 1, x = 0, y = 0, s = 0, t = 1, l = -2, v = 1)
easigma_bar <- function(j, x, y, s, t, l, v) {
    .Call(`_layeredBB_easigma_bar`, j, x, y, s, t, l, v)
}

#' Sigma (Equation 140 in ST329)
#'
#' This function evaluates the sigma function used to simulate Bessel Layers in the infinite sums
#'
#' @param j real value
#' @param x start value of Brownian bridge
#' @param y end value of Brownian bridge
#' @param s start value of Brownian bridge
#' @param t end value of Brownian bridge
#' @param l lower bound of Brownian bridge
#' @param v upper bound of Brownian bridge
#' 
#' @return real value: sigma evaluated at point j
#'
#' @examples
#' easigma(j = 1, x = 0, y = 0, s = 0, t = 1, l = -2, v = 1)
easigma <- function(j, x, y, s, t, l, v) {
    .Call(`_layeredBB_easigma`, j, x, y, s, t, l, v)
}

#' Phi_Bar (Equation 142 in ST329)
#'
#' This function evaluates the phi_bar function used to simulate Bessel Layers in the infinite sums
#'
#' @param j real value
#' @param x start value of Brownian bridge
#' @param y end value of Brownian bridge
#' @param s start value of Brownian bridge
#' @param t end value of Brownian bridge
#' @param l lower bound of Brownian bridge
#' @param v upper bound of Brownian bridge
#' 
#' @return real value: phi_bar evaluated at point j
#'
#' @examples
#' eaphi_bar(j = 1, x = 0, y = 0, s = 0, t = 1, l = -2, v = 1)
eaphi_bar <- function(j, x, y, s, t, l, v) {
    .Call(`_layeredBB_eaphi_bar`, j, x, y, s, t, l, v)
}

#' Phi (Equation 140 in ST329)
#'
#' This function evaluates the eaphi function used to simulate Bessel Layers in the infinite sums
#'
#' @param j real value
#' @param x start value of Brownian bridge
#' @param y end value of Brownian bridge
#' @param s start value of Brownian bridge
#' @param t end value of Brownian bridge
#' @param l lower bound of Brownian bridge
#' @param v upper bound of Brownian bridge
#' 
#' @return real value: eaphi evaluated at point j
#'
#' @examples
#' eaphi(j = 1, x = 0, y = 0, s = 0, t = 1, l = -2, v = 1)
eaphi <- function(j, x, y, s, t, l, v) {
    .Call(`_layeredBB_eaphi`, j, x, y, s, t, l, v)
}

#' Psi (Equation 155 and 162 in ST329)
#'
#' This function evaluates the eapsi function used to simulate Bessel Layers in the infinite sums
#'
#' @param j real value
#' @param xoy maximum of x and y where x and y are the start and end values of Brownian bridge
#' @param s start value of Brownian bridge
#' @param t end value of Brownian bridge
#' @param min minimum of Brownian bridge
#' @param v upper bound of Brownian bridge
#' 
#' @return real value: eapsi evaluated at point j
#'
#' @examples
#' eapsi(j = 1, xoy = 0, s = 0, t = 1, min = -2, v = 1)
eapsi <- function(j, xoy, s, t, min, v) {
    .Call(`_layeredBB_eapsi`, j, xoy, s, t, min, v)
}

#' Chi (Equation 156 and 163 in ST329)
#'
#' This function evaluates the chi function used to simulate Bessel Layers in the infinite sums
#'
#' @param j real value
#' @param xoy maximum of x and y where x and y are the start and end values of Brownian bridge
#' @param s start value of Brownian bridge
#' @param t end value of Brownian bridge
#' @param min minimum of Brownian bridge
#' @param v upper bound of Brownian bridge
#' 
#' @return real value: chi evaluated at point j
#'
#' @examples
#' eachi(j = 1, xoy = 0, s = 0, t = 1, min = -2, v = 1)
eachi <- function(j, xoy, s, t, min, v) {
    .Call(`_layeredBB_eachi`, j, xoy, s, t, min, v)
}

#' Gamma (Corollary 2 and Algorithm 26 in ST329) 
#'
#' This function evaluates the gamma function, S_{n}^{gamma},
#' used to simulate Bessel Layers in the infinite sums
#'
#' @param n integer value
#' @param x start value of Brownian bridge
#' @param y end value of Brownian bridge
#' @param s start value of Brownian bridge
#' @param t end value of Brownian bridge
#' @param l lower bound of Brownian bridge
#' @param v upper bound of Brownian bridge
#' 
#' @return real value: S_{n}^{gamma} evaluated at n
#'
#' @examples
#' eagamma(n = 1, x = 0, y = 0, s = 0, t = 1, l = -0.4, v = 0.8)
eagamma <- function(n, x, y, s, t, l, v) {
    .Call(`_layeredBB_eagamma`, n, x, y, s, t, l, v)
}

#' Delta_1 (Corollary 3 in ST329)
#'
#' This function evaluates the delta_1 function, S_{n}^{delta,1}
#' used to simulate Bessel Layers in the infinite sums
#'
#' @param n integer value
#' @param x start value of Brownian bridge
#' @param y end value of Brownian bridge
#' @param s start value of Brownian bridge
#' @param t end value of Brownian bridge
#' @param min minimum of Brownian bridge
#' @param v upper bound of Brownian bridge
#' 
#' @return real value: S_{n}^{delta,1} evaluated at n
#'
#' @examples
#' eadelta1(n = 1, x = 0, y = 0.8, s = 0, t = 1, min = -2, v = 2)
eadelta1 <- function(n, x, y, s, t, min, v) {
    .Call(`_layeredBB_eadelta1`, n, x, y, s, t, min, v)
}

#' Delta_2 (Corollary 4 in ST329)
#'
#' This function evaluates the delta_2function, S_{n}^{delta,2}
#' used to simulate Bessel Layers in the infinite sums
#'
#' @param n integer value
#' @param x start value of Brownian bridge
#' @param y end value of Brownian bridge
#' @param s start value of Brownian bridge
#' @param t end value of Brownian bridge
#' @param min minimum of Brownian bridge
#' @param v upper bound of Brownian bridge
#' 
#' @return real value: S_{n}^{delta,2} evaluated at n
#'
#' @examples
#' eadelta2(n = 1, x = -2, y = 0.8, s = 0, t = 1, min = -2, v = 2)
eadelta2 <- function(n, x, y, s, t, min, v) {
    .Call(`_layeredBB_eadelta2`, n, x, y, s, t, min, v)
}

#' Delta (Corollary 3 and 4 in ST329)
#'
#' This function evaluates the delta function, S_{n}^{delta}
#' used to simulate Bessel Layers in the infinite sums
#'
#' @param n integer value
#' @param x start value of Brownian bridge
#' @param y end value of Brownian bridge
#' @param s start value of Brownian bridge
#' @param t end value of Brownian bridge
#' @param min minimum of Brownian bridge
#' @param v upper bound of Brownian bridge
#' 
#' @return real value: S_{n}^{delta} evaluated at n
#'
#' @examples
#' # example where min(x,y) > min (delta_1 case)
#' eadelta(n = 1, x = 0,  y = 0.8, s = 0, t = 1, min = -2, v = 2)
#' eadelta1(n = 1, x = 0, y = 0.8, s = 0, t = 1, min = -2, v = 2)
#' # example where min(x,y) == min
#' eadelta(n = 1, x = -2, y = 0.8, s = 0, t = 1, min = -2, v = 2)
#' eadelta2(n = 1, x = -2, y = 0.8, s = 0, t = 1, min = -2, v = 2)
eadelta <- function(n, x, y, s, t, min, v) {
    .Call(`_layeredBB_eadelta`, n, x, y, s, t, min, v)
}

#' Calculate interval: [S^{gamma}_{2k+1}, S^{gamma}_{2k}] (Corollary 2 and Algorithm 26 in ST329)
#'
#' This function calculates the interval [S^{gamma}_{2k+1}, S^{gamma}_{2k}] for given k
#'
#' @param k integer value
#' @param x start value of Brownian bridge
#' @param y end value of Brownian bridge
#' @param s start value of Brownian bridge
#' @param t end value of Brownian bridge
#' @param l lower bound of Brownian bridge
#' @param v upper bound of Brownian bridge
#'
#' @return vector of two values, S^{gamma}_{2k+1} and S^{gamma}_{2k}
#'
#' @examples
#' eagamma_intervals(k = 1, x = 0, y = 0, s = 0, t = 1, l = -2, v = 1)
eagamma_intervals <- function(k, x, y, s, t, l, v) {
    .Call(`_layeredBB_eagamma_intervals`, k, x, y, s, t, l, v)
}

#' Calculate interval: [S^{delta,1}_{2k+1}, S^{delta,1}_{2k}] (Corollary 3 in ST329)
#'
#' This function calculates the interval [S^{delta,1}_{2k+1}, S^{delta, 1}_{2k}] 
#' (case where min(x,y) > min)
#'
#' @param k integer value
#' @param x start value of Brownian bridge
#' @param y end value of Brownian bridge
#' @param s start value of Brownian bridge
#' @param t end value of Brownian bridge
#' @param min minimum of Brownian bridge
#' @param v upper bound of Brownian bridge
#'
#' @return vector of two values, S^{delta,1}_{2k+1} and S^{delta,1}_{2k}
#'
#' @examples
#' eadelta1_intervals(k = 1, x = 0, y = 0, s = 0, t = 1, min = -2, v = 1)
eadelta1_intervals <- function(k, x, y, s, t, min, v) {
    .Call(`_layeredBB_eadelta1_intervals`, k, x, y, s, t, min, v)
}

#' Calculate interval: [S^{delta,2}_{2k+1}, S^{delta,2}_{2k}] (Corollary 4 in ST329)
#'
#' This function calculates the interval [S^{delta,2}_{2k+1}, S^{delta, 2}_{2k}] 
#' (case where min(x,y) == min)
#'
#' @param k integer value
#' @param x start value of Brownian bridge
#' @param y end value of Brownian bridge
#' @param s start value of Brownian bridge
#' @param t end value of Brownian bridge
#' @param min minimum of Brownian bridge
#' @param v upper bound of Brownian bridge
#'
#' @return vector of two values, S^{delta,2}_{2k+1} and S^{delta,2}_{2k}
#'
#' @examples
#' K = ceiling(sqrt((1)+(abs(1-(-2))*abs(1-(-2))))/(2*abs(1-(-2))))
#' eadelta2_intervals(k = K, x = -2, y = 0, s = 0, t = 0, min = -2, v = 1)
eadelta2_intervals <- function(k, x, y, s, t, min, v) {
    .Call(`_layeredBB_eadelta2_intervals`, k, x, y, s, t, min, v)
}

#' Calculate interval: [S^{delta}_{2k+1}, S^{delta}_{2k}] (Algorithm 28 in ST329)
#'
#' This function calculates the interval [S^{delta}_{2k+1}, S^{delta}_{2k}] 
#' (case where min(x,y) > min or where min(x,y) == min)
#'
#' @param k integer value
#' @param x start value of Brownian bridge
#' @param y end value of Brownian bridge
#' @param s start value of Brownian bridge
#' @param t end value of Brownian bridge
#' @param min minimum of Brownian bridge
#' @param v upper bound of Brownian bridge
#'
#' @return vector of two values, S^{delta}_{2k+1} and S^{delta}_{2k}
#'
#' @examples
#' # case where min(x,y ) > min
#' eadelta1_intervals(k = 1, x = 0, y = 0, s = 0, t = 1, min = -2, v = 1)
#'
#' # case where min(x,y) == min
#' K = ceiling(sqrt((1)+(abs(1-(-2))*abs(1-(-2))))/(2*abs(1-(-2))))
#' eadelta_intervals(k = K, x = -2, y = 0, s = 0, t = 0, min = -2, v = 1)
eadelta_intervals <- function(k, x, y, s, t, min, v) {
    .Call(`_layeredBB_eadelta_intervals`, k, x, y, s, t, min, v)
}

#' Find product of a vector
#'
#' This function product of the elements of a numerical vector
#'
#' @param vect numerical vector
#'
#' @return product of element in vector given
#'
#' @examples
#' # returns 120
#' product_vector(c(1,2,3,4,5)) 
product_vector <- function(vect) {
    .Call(`_layeredBB_product_vector`, vect)
}

#' Gamma coin flipper (Algorithm 26 in ST329)
#'
#' Flips 'Gamma coin'; uses the Cauchy sequence S^{gamma}_{k} to 
#' determine whether or not the Brownian bridge starting at x, ending at y, between [s,t]
#' remains in interval [l,v]
#'
#' @param u simulated value from random U[0,1] 
#' @param k integer value starting index for calculating the intervals
#' @param x start value of Brownian bridge
#' @param y end value of Brownian bridge
#' @param s start value of Brownian bridge
#' @param t end value of Brownian bridge
#' @param l lower bound of Brownian bridge
#' @param v upper bound of Brownian bridge
#'
#' @examples
#' gamma_coin(u = runif(1, 0, 1),
#'            k = 0,
#'            x = 0,
#'            y = 0,
#'            s = 0,
#'            t = 1,
#'            l = -0.5,
#'            v = 0.5)
#'
#' @return boolean value: if T, accept probability that Brownian bridge remains 
#'         in [l,v], otherwise reject
gamma_coin <- function(u, k, x, y, s, t, l, v) {
    .Call(`_layeredBB_gamma_coin`, u, k, x, y, s, t, l, v)
}

#' Gamma coin flipper for intervals
#'
#' Flips 'Gamma coin' for intervals; takes the product of the Cauchy sequence S^{gamma}_{k} to
#' determine whether or not the Brownian bridge remains in the interval [l,v]
#'
#' @param u simulated value from random U[0,1] 
#' @param k integer value starting index for calculating the intervals
#' @param X vector of values of Brownian bridge
#' @param times vector of times
#' @param l lower bound of Brownian bridge
#' @param v upper bound of Brownian bridge
#'
#' @examples
#' # setting up Brownian bridge variable
#' brownian_bridge <- matrix(c(0, 0, -0.2, 0.4, 0.3, 0.5, 1, 1),
#'                           ncol = 4, nrow = 2)
#' 
#' # flip delta coin whether or not Brownian bridge remains in [-0.5, 1.5]
#' gamma_coin_intervals(u = runif(1, 0, 1),
#'                      k = 1,
#'                      X = brownian_bridge[1,],
#'                      times = brownian_bridge[2,],
#'                      l = -0.5,
#'                      v = 1.5)
#'
#' @return boolean value: if T, accept probability that Brownian bridge remains 
#'         in [l,v], otherwise reject
gamma_coin_intervals <- function(u, k, X, times, l, v) {
    .Call(`_layeredBB_gamma_coin_intervals`, u, k, X, times, l, v)
}

#' Delta coin flipper (Algorithm 28 in ST329)
#'
#' Flips 'Delta coin'; uses the Cauchy sequence S^{delta}_{k} to 
#' determine whether or not the Brownian bridge with minimum, min, 
#' starting at x, ending at y, between [s,t] remains in interval [l,v]
#'
#' @param u simulated value from random U[0,1] 
#' @param k integer value starting index for calculating the intervals
#' @param x start value of Brownian bridge
#' @param y end value of Brownian bridge
#' @param s start value of Brownian bridge
#' @param t end value of Brownian bridge
#' @param min minimum of Brownian bridge
#' @param v upper bound of Brownian bridge
#'
#' @examples 
#' delta_coin(u = runif(1, 0, 1),
#'            k = 0,
#'            x = 0.1,
#'            y = 0.4,
#'            s = 0,
#'            t = 1,
#'            min = -0.2,
#'            v = 1.5)
#'                      
#' @return boolean value: if T, accept probability that Brownian bridge with 
#'         minimum, min, remains in [l,v], otherwise reject
delta_coin <- function(u, k, x, y, s, t, min, v) {
    .Call(`_layeredBB_delta_coin`, u, k, x, y, s, t, min, v)
}

#' Delta coin flipper for intervals (used for Algorithm 33 in ST329)
#'
#' Flips 'Delta coin' for intervals; takes the product of the Cauchy sequence S^{delta}_{k} to 
#' determine whether or not the Brownian bridge with minimum, min, remains in the interval [l,v]
#'
#' @param u simulated value from random U[0,1] 
#' @param k integer value starting index for calculating the intervals
#' @param X vector of values of Brownian bridge
#' @param times vector of times
#' @param min minimum of Brownian bridge
#' @param v upper bound of Brownian bridge
#'
#' @examples
#' # setting up Brownian bridge variable
#' brownian_bridge <- matrix(c(0, 0, -0.2, 0.4, 0.3, 0.5, 1, 1),
#'                           ncol = 4, nrow = 2)
#' 
#' # flip delta coin whether or not Brownian bridge remains in [-0.2, 1.5]
#' d <- abs(1.5 - -0.2)
#' k <- ceiling(sqrt(1 + d^2)/(2*d))
#' delta_coin_intervals(u = runif(1, 0, 1),
#'                      k = k,
#'                      X = brownian_bridge[1,],
#'                      times = brownian_bridge[2,],
#'                      min = -0.2,
#'                      v = 1.5)
#'
#' @return boolean value: if T, accept probability that Brownian bridge with 
#'         minimum, min, remains in [min,v], otherwise reject
delta_coin_intervals <- function(u, k, X, times, min, v) {
    .Call(`_layeredBB_delta_coin_intervals`, u, k, X, times, min, v)
}

#' Inverse Gaussian Sampler
#'
#' This function returns 1 sample from an Inverse Gaussian 
#' distribution with mean mu and shape lambda
#'
#' @param mu mean value
#' @param lambda shape value
#' 
#' @return real value: simulated point from Inverse Gaussian 
#'         distribution with mean mu and shape lambda
#'
#' @examples
#' curve(statmod::dinvgauss(x, mean = 1, shape = 1), 0, 4)
#' samples <- sapply(1:10000, function(i) inv_gauss_sampler(mu = 1, lambda = 1))
#' lines(density(x = samples, adjust = 0.5), col = 'blue')
inv_gauss_sampler <- function(mu, lambda) {
    .Call(`_layeredBB_inv_gauss_sampler`, mu, lambda)
}

#' Layered Brownian Bridge sampler (Algorithm 33 in ST329)
#'
#' Simulation of a layered Brownian Bridge given a Bessel layer at user specified times
#'
#' @param x start value of Brownian bridge
#' @param y end value of Brownian bridge
#' @param s start time of Brownian bridge
#' @param t end time of Brownian bridge
#' @param a vector/sequence of numbers
#' @param l integer number denoting Bessel layer, i.e. Brownian bridge 
#'        is contained in [min(x,y)-a[l], max(x,y)+a[l]]
#' @param times vector of real numbers to simulate Bessel bridge
#' 
#' @return A list with the following components
#' \describe{
#'   \item{full_path}{Matrix of the simulated layered Brownian bridge path at 
#'                    all included time points, i.e. s, t and times. The times
#'                    are sorted and duplicates are removed. The first row
#'                    are the points of the Brownian bridge (named 'X') 
#'                    second row are corresponding times (named 'time')}
#'   \item{simulated_path}{Matrix of the simulated layered Brownian bridge path 
#'                         only at the specified times passed into the function, 
#'                         i.e. the times vector. The times are not sorted and 
#'                         duplicates are not removed. The first row 
#'                         are the points of the layered Brownian bridge (named 
#'                         'X') second row are corresponding times (named 
#'                         'times')}
#'   \item{remove_m_path}{Matrix of the simulated layered Brownian bridge path 
#'                        only at all included times points excluding tau. These 
#'                        times are sorted and duplicates are removed. The first
#'                        row are the points of the layered Brownian bridge 
#'                        (named 'X') second row are corresponding times 
#'                        (named 'time')}
#' }
#'
#' @examples
#' # simulate Bessel layer
#' bes_layer <- bessel_layer_simulation(x = 0,
#'                                      y = 0,
#'                                      s = 0,
#'                                      t = 1,
#'                                      mult = 0.2)
#' # simulate layered Brownian bridge
#' # notice full_path has all included times and are sorted and have no duplicates
#' # simulated_path only returns points that are passed into times
#' # remove_m does not include the simulated minimum or maximum point
#' layered_brownian_bridge(x = 0,
#'                         y = 0,
#'                         s = 0,
#'                         t = 1,
#'                         bessel_layer = bes_layer,
#'                         times = c(0.2, 0.4, 0.6, 0.8))
#'
#' # note that simulated_path does not remove duplicates passed into times
#' layered_brownian_bridge(x = 0,
#'                         y = 0,
#'                         s = 0,
#'                         t = 1,
#'                         bessel_layer = bes_layer,
#'                         times = c(0.2, 0.4, 0.6, 0.8, 0.4, 0.6))
#'
#' # another example
#' start <- runif(1, -1, 1)
#' end <- runif(1, -1, 1)
#' bes_layer <- bessel_layer_simulation(x = start, y = end, s = 0, t = 1, mult = 0.2)
#' path <- layered_brownian_bridge(x = start,
#'                                 y = end,
#'                                 s = 0,
#'                                 t = 1,
#'                                 bessel_layer = bes_layer,
#'                                 times = seq(0, 1, 0.01))$full_path
#' plot(x = path['time',], y = path['X',], pch = 20, xlab = 'Time', ylab = 'X',
#'      ylim = c(bes_layer$L, bes_layer$U))
#' lines(x = path['time',], y = path['X',])
#' abline(h=c(bes_layer$L, bes_layer$U), col = 'red')
#' abline(h=c(bes_layer$l, bes_layer$u), col = 'red', lty = 2)
#' 
#' # compare the simulated distribution of simulated points to the
#' # theoretical distribution of simulated points
#' # for large Bessel layers, it should look like a unconditional Brownian bridge
#' x <- 0.53
#' y <- 4.32
#' s <- 0.53
#' t <- 2.91
#' q <- 1.72
#' replicates <- 10000
#' paths <- list()
#' large_bessel_layer <- bessel_layer_simulation(x = x,
#'                                               y = y,
#'                                               s = s,
#'                                               t = t,
#'                                               mult = 100)
#' # repeatedly simulate Brownian bridge 
#' for (i in 1:replicates) {
#'   paths[[i]] <- layered_brownian_bridge(x = x,
#'                                         y = y,
#'                                         s = s,
#'                                         t = t,
#'                                         bessel_layer = large_bessel_layer,
#'                                         times = seq(s, t, 0.01))
#' }
#' # select the points at the specified time q
#' index <- which(seq(s, t, 0.01)==q)
#' simulated_points <- sapply(1:replicates, function(i) paths[[i]]$full_path['X', index])
#' # calculate the theoretical mean and standard deviation of the simulated points at time q
#' theoretical_mean <- x + (q-s)*(y-x)/(t-s)
#' theoretical_sd <- sqrt((t-q)*(q-s)/(t-s))
#' # plot distribution of the simulated points and the theoretical distribution
#' plot(density(simulated_points))
#' curve(dnorm(x, theoretical_mean, theoretical_sd), add = T, col = 'red')
layered_brownian_bridge <- function(x, y, s, t, bessel_layer, times) {
    .Call(`_layeredBB_layered_brownian_bridge`, x, y, s, t, bessel_layer, times)
}

#' Multi-dimensional Layered Brownian Bridge sampler
#'
#' Simulation of a multi-dimensional layered Brownian Bridge given Bessel layers
#' at user-specified times
#'
#' @param dim dimension of Brownian bridge
#' @param x start value of Brownian bridge
#' @param y end value of Brownian bridge
#' @param s start time of Brownian bridge
#' @param t end time of Brownian bridge
#' @param bessel_layers a list of length dim where list[i] is the Bessel layer for component i
#' @param times vector of real numbers to simulate Bessel bridge
#' 
#' @return A list with the following components
#' \describe{
#'   \item{full_path}{Matrix of the simulated layered Brownian bridge path at all 
#'                    included time points, i.e. s, t and times. The times
#'                    are sorted and duplicates are removed. The first dim rows
#'                    are the points of the Brownian bridge in each component,
#'                    last row gives the corresponding times}
#'   \item{simulated_path}{Matrix of the simulated layered Brownian bridge path 
#'                         only at the specified times passed into the function, 
#'                         i.e. the times vector. The times are not sorted and
#'                         duplicates are not removed. The first dim rows
#'                         are the points of the layered Brownian bridge in 
#'                         each component, last row gives the corresponding times}
#' }
#'
#' @examples
#' # simulate Bessel layer for two-dimensional Brownian bridge starting 
#' # and ending at (0,0) in time [0,1]
#' bes_layers <- multi_bessel_layer_simulation(dim = 2,
#'                                             x = c(0, 0),
#'                                             y = c(0, 0),
#'                                             s = 0,
#'                                             t = 1,
#'                                             mult = 0.5)
#' # simulate two-dimensional Brownian bridge starting 
#' # and ending at (0,0) in time [0,1]
#' multi_layered_brownian_bridge(dim = 2,
#'                               x = c(0,0),
#'                               y = c(0,0),
#'                               s = 0,
#'                               t = 1,
#'                               bessel_layers = bes_layers,
#'                               times = c(0.2, 0.4, 0.6, 0.8))
#'
#' # note that simulated_path does not remove duplicates passed into times
#' multi_layered_brownian_bridge(dim = 2,
#'                               x = c(0,0),
#'                               y = c(0,0),
#'                               s = 0,
#'                               t = 1,
#'                               bessel_layers = bes_layers,
#'                               times = c(0.2, 0.4, 0.6, 0.8, 0.4, 0.6))
multi_layered_brownian_bridge <- function(dim, x, y, s, t, bessel_layers, times) {
    .Call(`_layeredBB_multi_layered_brownian_bridge`, dim, x, y, s, t, bessel_layers, times)
}