|
11 | 11 | #' |
12 | 12 | #' @return a function implementing the derivative as a finite-difference approximation. |
13 | 13 | #' This has a second argument, `.h`, that allow the finite-difference to be set when evaluating |
14 | | -#' the function. The default values are set for reasonable numerical precision. |
| 14 | +#' the function. The default values are set for reasonable numerical precision |
| 15 | +#' with the pattern-book functions. |
15 | 16 | #' |
16 | 17 | #' @details |
17 | 18 | #' Uses a simple finite-difference scheme to evaluate the derivative. The function created |
@@ -64,11 +65,15 @@ numD <- function(tilde, ..., .h=NULL) { |
64 | 65 | the_h <- ifelse(is.null(.h), 0.000001, .h) |
65 | 66 | res = make_dfdx( f, dvars[1], .h = the_h) %>% |
66 | 67 | bind_params(formals(f)) |
67 | | - }else if (length(dvars==2) && dvars[1]==dvars[2]) { |
| 68 | + } else if (length(dvars)==2 && dvars[1]==dvars[2]) { |
68 | 69 | # Second unmixed partial |
69 | 70 | the_h <- ifelse(is.null(.h), 0.0001, .h) |
70 | 71 | res = make_d2fd2x( f, dvars[1], .h = the_h) %>% |
71 | 72 | bind_params(formals(f)) |
| 73 | + } else if (length(dvars)==3 && length(unique(dvars)) == 1) { |
| 74 | + the_h = ifelse(is.null(.h), 0.01, .h) |
| 75 | + res = make_d3fd3x( f, dvars[1], .h = the_h) %>% |
| 76 | + bind_params(formals(f)) |
72 | 77 | } else if (length(dvars)==2) { |
73 | 78 | # mixed partial |
74 | 79 | the_h <- ifelse(is.null(.h), 0.001, .h) |
@@ -143,3 +148,27 @@ make_d2fd2x <- function(f, .wrt, .h = 0.000001) { |
143 | 148 |
|
144 | 149 | conventional_argument_order(dfun, ".h") |
145 | 150 | } |
| 151 | + |
| 152 | +# ========== |
| 153 | +# |
| 154 | +# @note Helper function for third-order deriv in one variable. |
| 155 | + |
| 156 | +make_d3fd3x <- function(f, .wrt, .h = 0.0001) { |
| 157 | + far_right_args <- far_left_args <- right_args <- left_args <- args <- names(formals(f)) |
| 158 | + args <- paste0(args, collapse=", ") |
| 159 | + right_args[right_args==.wrt] <- glue::glue("{.wrt} + .h") |
| 160 | + right_args <- paste0(right_args, collapse=", ") |
| 161 | + left_args[left_args==.wrt] <- glue::glue("{.wrt} - .h") |
| 162 | + left_args <- paste0(left_args, collapse=", ") |
| 163 | + far_right_args[far_right_args==.wrt] <- glue::glue("{.wrt} + 2*.h") |
| 164 | + far_right_args <- paste0(far_right_args, collapse=", ") |
| 165 | + far_left_args[far_left_args==.wrt] <- glue::glue("{.wrt} - 2*.h") |
| 166 | + far_left_args <- paste0(far_left_args, collapse=", ") |
| 167 | + command <- glue::glue( |
| 168 | + "function({.wrt}){{(f({far_right_args}) - 2*(f({right_args}) - f({left_args})) - f({far_left_args}))/(2*.h^3)}}") |
| 169 | + dfun <- eval(parse(text=command)) |
| 170 | + formals(dfun) <- c(formals(f), list(.h=.h)) |
| 171 | + |
| 172 | + conventional_argument_order(dfun, ".h") |
| 173 | +} |
| 174 | + |
0 commit comments