R, Monte Carlo y problemas empresariales, parte 2

Por paradójico que parezca, pero todavía con bastante frecuencia en la empresa hay tareas que son diferentes a construir otra cuenta personal, otro monitoreo u otro flujo de trabajo. Si piensa un poco y no se aferra a codificar o buscar software especializado, puede escribir una solución compacta, muy elegante y rápida utilizando el método Monte Carlo.







Las tareas empresariales son lo suficientemente compactas como para iterar y no requieren 100 lugares decimales. No estamos lanzando cohetes o reactores y construyendo una teoría científica de todo.







Consideremos más adelante el ejemplo de una de las tareas no estándar.







Es una continuación de una serie de publicaciones anteriores .







Formulación del problema



El desafío tiene sus raíces en IoT para la agricultura. Es decir, la colocación de sensores en un campo de papa con un patrón de riego circular. Pidamos a Google algunas imágenes: "La respuesta al misterio de los campos redondos: ¡todo es más interesante de lo que crees!" ... Es necesario proporcionar las características necesarias de la cobertura de la red en malla, teniendo en cuenta las distancias permisibles entre vecinos. Al mismo tiempo, es necesario optimizar el presupuesto y emitir coordenadas gps para la colocación de sensores y formar el esquema de bypass más corto.







Decisión



Conectamos paquetes



Paquetes
library(tidyverse)
library(magrittr)
library(scales)
library(data.table)
library(tictoc)
library(stringi)
library(arrangements)
library(ggforce)
      
      





Descomposición



.







  1. N



    .
  2. .
  3. , N



    1. .
  4. .


. ? -.







-



.







drawDisk <- function(df) {
  #      
  #    ,      0
  for(col in c("force_x", "force_y")){
    if (!(col %in% names(df))) df[, col] <- 0
  }

  ggplot(data = df, aes(x = x, y = y)) +
    ggforce::geom_circle(aes(x0 = 0, y0 = 0, r = 1, colour = "red"), 
                         inherit.aes = FALSE) +
    geom_point(size = 2, fill = "black", shape = 21) +
    geom_text(aes(label = idx), hjust = 0.5, vjust = -1) +
    #   
    geom_segment(aes(xend = x + force_x / 10, yend = y + force_y / 10), 
                 colour = "red", 
                 arrow = arrow(length = unit(0.2,"cm")), size = 0.6) + 
    xlim(-1.5, 1.5) +
    ylim(-1.5, 1.5) +
    coord_fixed() +
    labs(x = " X", y = " Y") +
    theme_bw()
}
      
      







, -. , . ( ) .







— .







, N



. , . ( ). . , . , , «». .







#  -    .
# ,     ,   
charges_dt <- tibble(idx = 1:13) %>%
  mutate(angle = runif(n(), min = 0, max = 2*pi), 
         r = runif(n(), min = 0, max = 1),
         x = r * cos(angle), y = r * sin(angle)) %>%
  select(idx, x, y) %>%
  setDT(key = "idx")

#          
keepers_dt <- max(charges_dt$idx) %>% 
  {tibble(idx = (. + 1):(. + 40))} %>%
  mutate(angle = (idx - 1) * (2 * pi / n()),
         x = 1.3 * cos(angle), y = 1.3 * sin(angle)) %>%
  select(idx, x, y) %>%
  setDT(key = "idx")
      
      





, .







full_dt <- rbindlist(list(charges_dt, keepers_dt))
drawDisk(full_dt)
      
      







, c ( ).

:







  • (- );
  • .






s=at2/2=(F/m)t2/2

















Δs=FΔt







, .







max_force <- 10

tic("Balancing charges")
#   !
#        
#      0 
#       
while (max_force > 0.05 & nrow(charges_dt[x^2 + y^2 > 1.2]) == 0) {
  #       
  full_dt <- rbindlist(list(charges_dt, keepers_dt), fill = TRUE)

  ff <- function(x0, y0){
    #   -- 1/r2,  ;
    #    , sqrt(r2) -- 
    #       

    dx = full_dt$x - x0
    dy = full_dt$y - y0
    r2 = dx^2 + dy^2
    # na.rm  NaN  ..
    list(sum(-dx * r2^(-1.5), na.rm = TRUE),
         sum(-dy * r2^(-1.5), na.rm = TRUE))
  }  

  #   ,    " " 
  charges_dt[, c("force_x", "force_y") := ff(x0 = x, y0 = y), by = idx]
  #   ,   
  max_force <- charges_dt[, max(sqrt(force_x^2 + force_y^2), na.rm = TRUE)]
  force_scale = if_else(max_force > 1, 1 / max_force / 1e2, 1/ max_force / 5e2)

  #   
  charges_dt %>%
    .[, `:=`(x = x + force_x * force_scale, 
             y = y + force_y * force_scale)]
}
toc()

full_dt <- rbindlist(list(charges_dt, keepers_dt), fill = TRUE)
      
      







-. :







  • ;
  • , , ( );
  • .


optimizePath <- function(dt) {
  #       
  # 1.      ,      
  dt[, r0 := sqrt(x^2+y^2)] %>%
    setorder(-r0)
  n1 <- dt[1, idx]

  #       
  #      n1,     
  points_in <- dt[idx != n1, idx]

  #         
  # (  ,   )
  augment_tbl <- dt %>%
    mutate_at("idx", `*`, -1) %>%
    mutate(r0 = sqrt(x^2 + y^2)) %>%
    mutate_at(vars(x, y), ~(.x/r0)) %>%
    bind_rows(dt) %>%
    select(idx, x, y)

  #      
  ll_tbl <- unique(augment_tbl$idx) %>%
    tidyr::expand_grid(l = ., r = .) %>%
    filter(l != r, (l > 0 | r > 0)) %>%
    #    
    rowwise() %>%
    # mutate(m = list(sort(c(l, r))))
    mutate(edge_id = stri_c(sort(c(l, r)), collapse = "=")) %>%
    ungroup() %>%
    distinct(edge_id, .keep_all = TRUE) %>%
    #    
    left_join(select(augment_tbl, idx, l_x = x, l_y = y), by = c("l" = "idx")) %>%
    #    
    left_join(select(augment_tbl, idx, r_x = x, r_y = y), by = c("r" = "idx")) %>%
    mutate(s = sqrt((l_x - r_x)^2 + (l_y - r_y)^2)) %>%
    arrange(l, r)

  points_seq <- arrangements::permutations(v = points_in, replace = FALSE, 
                                       layout = "column", nsample = 5000)
  #        .     
  routes_lst <- points_seq %>% 
    rbind(-n1, n1, ., -tail(., 1)) %>%
    as.data.frame() %>% as.list()

  #       
  routes_dt <- data.table(route_id = seq_along(routes_lst), route = routes_lst) %>%
    .[, .(from = unlist(route)), by = route_id] %>%
    .[, to := shift(from, type = "lead"), by = route_id] %>%
    #    
    na.omit() %>%
    #    
    .[, edge_id := stri_c(sort(unlist(.BY)), collapse = "="), by = .(from, to)] %>%
    .[, .(route_id, edge_id)] %>%
    #       
    .[as.data.table(ll_tbl), s := i.s, on = "edge_id"]

  #   ,  
  best_routes <- routes_dt[, .(len = sum(s)), by = route_id] %>%
    setorder(len) %>%
    head(10) %T>%
    print()

  #  -10  
  best_routes %>%
    select(route_id) %>%
    mutate(idx = routes_lst[route_id]) %>%
    tidyr::unnest(idx) %>%
    left_join(augment_tbl) %>%
    tidyr::nest(data = -route_id) %>%
    left_join(best_routes)
}
      
      











    route_id      len
 1:     2070 8.332854
 2:     2167 8.377680
 3:     4067 8.384417
 4:     3614 8.418678
 5:     5000 8.471521
 6:     4495 8.542041
 7:     2233 8.598278
 8:     4430 8.609391
 9:     2915 8.616048
10:     3380 8.695547
      
      











tic("Route optimization")
best_tbl <- optimizePath(charges_dt)
toc()

best_route_tbl <- best_tbl$data[[1]]
full_dt <- rbindlist(list(best_route_tbl, keepers_dt), fill = TRUE)
gp <- drawDisk(full_dt) +
  #   
  geom_path(arrow = arrow(type = "closed"), data = best_route_tbl)

gp
      
      





Ruta de circunvalación









, . GPS, GPS . .









, . « », . tidyverse



$10^3$-$10^4$ . . , , C++. , .







  1. data.table



    . Base R .
  2. , . .
  3. data.frame



    . . tidyverse



    .
  4. , , .
  5. Monte Carlo es un muy buen enfoque. Una aplicación inicial rápida puede dar un resultado útil, así como buscar la solución al problema y, posiblemente, encontrar algunas simplificaciones y análisis.
  6. No dude en utilizar métodos de analogía. Pueden hacer posible construir un modelo simplificado del problema original, que es computacionalmente mucho más simple que el original y se puede transferir fácilmente a Monte Carlo.


Publicación anterior - "Niños, ruso y R" .








All Articles