diff --git a/R/si_model.R b/R/si_model.R
index e0a5c51..b12cd73 100644
--- a/R/si_model.R
+++ b/R/si_model.R
@@ -25,13 +25,13 @@
 #'   n = destination_all, d = distance_euclidean)
 #' head(od_output)
 #' plot(od$distance_euclidean, od_output$interaction)
-#' od_pconst = si_calculate(od, fun = fun, beta = 0.3, O = origin_all, n = destination_all,
-#'   d = distance_euclidean, constraint_p = origin_all)
+#' od_pconst = si_calculate(od, fun = fun, beta = 0.3, O = origin_all,
+#'   n = destination_all, d = distance_euclidean, constraint_p = origin_all)
 #' plot(od_pconst$distance_euclidean, od_pconst$interaction)
 #' plot(od_pconst["interaction"], logz = TRUE)
 si_calculate = function(od, fun, constraint_p, ...) {
   dots = rlang::enquos(...)
-  od = mutate(od, interaction = fun(!!!dots))
+  od = dplyr::mutate(od, interaction = fun(!!!dots))
   if (!missing(constraint_p)) {
     od_grouped = dplyr::group_by(od, .data$O)
     od_grouped = dplyr::mutate(
diff --git a/man/si_calculate.Rd b/man/si_calculate.Rd
index c1fd783..b917b68 100644
--- a/man/si_calculate.Rd
+++ b/man/si_calculate.Rd
@@ -26,16 +26,17 @@ and user-specified function
 }
 \examples{
 od = si_to_od(si_zones, si_zones, max_dist = 4000)
-fun_dd = function(od, d = "distance_euclidean", beta = 0.3) exp(-beta * od[[d]] / 1000)
-od_dd = si_calculate(od, fun = fun_dd)
+fun_dd = function(d = "distance_euclidean", beta = 0.3) exp(-beta * d / 1000)
+fun_dd(d = (1:5)*1000)
+od_dd = si_calculate(od, fun = fun_dd, d = distance_euclidean)
 plot(od$distance_euclidean, od_dd$interaction)
-fun = function(od, O, n, d, beta) od[[O]] * od[[n]] * exp(-beta * od[[d]]/1000)
-od_output = si_calculate(od, fun = fun, beta = 0.3, O = "origin_all", 
-  n = "destination_all", d = "distance_euclidean")
+fun = function(O, n, d, beta) O * n * exp(-beta * d / 1000)
+od_output = si_calculate(od, fun = fun, beta = 0.3, O = origin_all, 
+  n = destination_all, d = distance_euclidean)
 head(od_output)
 plot(od$distance_euclidean, od_output$interaction)
-od_pconst = si_calculate(od, fun = fun, beta = 0.3, O = "origin_all", n = "destination_all",
-  d = "distance_euclidean", constraint_p = origin_all)
+od_pconst = si_calculate(od, fun = fun, beta = 0.3, O = origin_all,
+  n = destination_all, d = distance_euclidean, constraint_p = origin_all)
 plot(od_pconst$distance_euclidean, od_pconst$interaction)
 plot(od_pconst["interaction"], logz = TRUE)
 }
diff --git a/vignettes/si.Rmd b/vignettes/si.Rmd
index aec4b99..3bfda1e 100644
--- a/vignettes/si.Rmd
+++ b/vignettes/si.Rmd
@@ -119,10 +119,13 @@ Note that the output has duplicate columns: `si_to_od()` joins data from the ori
 A simplistic SIM can be created just based on the distance between points:
 
 ```{r unconstrained1}
-si_power = function(od, beta) {
-  (od[["distance_euclidean"]] / 1000)^beta
-}
-od_calculated = si_calculate(od_sim, fun = si_power, beta = -0.8)
+si_power = function(d, beta) (d / 1000)^beta
+od_calculated = si_calculate(
+  od_sim,
+  fun = si_power,
+  d = distance_euclidean,
+  beta = -0.8
+  )
 plot(od_calculated["interaction"], logz = TRUE)
 ```
 
@@ -178,6 +181,7 @@ od_originating = od_joined %>%
 od_constrained_p = si_calculate(
   od_originating,
   fun = si_power,
+  d = distance_euclidean,
   beta = -0.8,
   constraint_p = originating_per_zone
   )