Merge pull request #157 from aadler/master

Reverse direction of hin/hinjac in wrapper functions and update various examples.
astamm · Jun 15, 2024 · 68763bc · 68763bc
2 parents e56b287 + fe63b95
commit 68763bc
Show file tree

Hide file tree

Showing 26 changed files with 1,340 additions and 977 deletions.
diff --git a/DESCRIPTION b/DESCRIPTION
@@ -1,7 +1,7 @@
 Package: nloptr
 Type: Package
 Title: R Interface to NLopt
-Version: 2.0.4
+Version: 2.1.0
 Authors@R: c(person("Jelmer", "Ypma", role = "aut",
                     email = "[email protected]"),
              person(c("Steven", "G."), "Johnson", role = "aut",

diff --git a/NEWS.md b/NEWS.md
@@ -1,14 +1,27 @@
-# nloptr 2.0.4
-
-* Updated roxygen version;
-* Updated maintainer email;
+# nloptr 2.1.0
+This release deprecates the default behavior of the inequality equations in any
+wrapper function which uses them. Currently, they are calibrated to be >= 0.
+This version allows for the equations to be consistent with the main `nloptr`
+function, which requires <= 0. In a future release, the default behavior will
+switch to assuming the calibration is <= 0, and eventually, the >= 0 behavior
+will be removed. It also includes a large number of safety and efficiency
+changes, and an expansion of the unit tests to 100% coverage for all files but
+one. The major changes include:
+
+* Reversed the direction of the inequality equations `hin` and `hinjac` in the
+wrapper functions which use them, bringing them into compliance with the main
+`nloptr` call. This addresses
+[Issue #148](https://github.com/astamm/nloptr/issues/148).
+* Cleaned the Hock-Schittkowski problem no. 100, Hartmann 6-dimensional, and
+Powell exponential examples. This addresses
+[Issue #152](https://github.com/astamm/nloptr/issues/152) and
+[Issue #156](https://github.com/astamm/nloptr/issues/156).
+* Updated roxygen version.
+* Updated maintainer email.
 * Deal with NA returns from detectCords (contributed by @jeroen in PR #150);
-* Setup rhub v2 checks;
-* Update cmake installation instructions on Mac with brew (#146);
-* Allow use of equality constraints with COBYLA (#135);
-
-# nloptr 2.0.3.9100
-
+* Setup rhub v2 checks.
+* Update cmake installation instructions on Mac with brew (#146).
+* Allow use of equality constraints with COBYLA (#135).
 * Replaced the unit testing framework of `testthat` with `tinytest` (See 
 [Issue #136](https://github.com/astamm/nloptr/issues/136)).
 * Brought coverage of `is.nloptr` to 100%. The only file not completely covered
@@ -17,18 +30,6 @@ trapped by tests in R before the call gets to C.
 * Linted package for code correctness and consistency.
 * Updated vignette, DESCRIPTION, and NEWS.
 * Updated package website to use bootstrap 5.
-
-# nloptr 2.0.3.9000
-
-This is a patch version update which should make the code safer, more efficient,
-and easier to follow. Please see commit logs for
-[#128](https://github.com/astamm/nloptr/pull/128),
-[#129](https://github.com/astamm/nloptr/pull/129),
-[#131](https://github.com/astamm/nloptr/pull/131), 
-[#132](https://github.com/astamm/nloptr/pull/132),
-and [#133](https://github.com/astamm/nloptr/pull/133) for the full
-description of the changes which include:
-
 * Expanded unit tests: coverage now over 97% with no file below 90%
 * Removed forcing `C++11`
 * Added safety checks to C code
@@ -38,6 +39,8 @@ description of the changes which include:
 * Updated Github actions
 * Some bugfixes (e.g. in `isres` or the warning in `nl.grad`.)
 
+Please see the commit logs for more detailed descriptions of the changes.
+
 # nloptr 2.0.3
 
 * Improved compatibility on RHEL/CentOS by first searching for a `cmake3` binary

diff --git a/R/auglag.R b/R/auglag.R
@@ -7,10 +7,15 @@
 #
 # Wrapper to solve optimization problem using Augmented Lagrangian.
 #
-# Changelog:
-#   2017-09-26: Fixed bug, BOBYQA is allowed as local solver
-#         (thanks to Leo Belzile).
-#   2023-02-08: Tweaks for efficiency and readability (Avraham Adler)
+# CHANGELOG
+#
+# 2017-09-26: Fixed bug, BOBYQA is allowed as local solver
+#   (thanks to Leo Belzile).
+# 2023-02-08: Tweaks for efficiency and readability (Avraham Adler)
+# 2024-06-04: Switched desired direction of the hin/hinjac inequalities, leaving
+#   the old behavior as the default for now. Also corrected Powell example.
+#   (Avraham Adler)
+#
 
 #' Augmented Lagrangian Algorithm
 #'
@@ -54,6 +59,10 @@
 #' the local solver?; not possible at the moment.
 #' @param nl.info logical; shall the original NLopt info been shown.
 #' @param control list of options, see \code{nl.opts} for help.
+#' @param deprecatedBehavior logical; if \code{TRUE} (default for now), the old
+#' behavior of the Jacobian function is used, where the equality is \eqn{\ge 0}
+#' instead of \eqn{\le 0}. This will be reversed in a future release and
+#' eventually removed.
 #' @param ... additional arguments passed to the function.
 #'
 #' @return List with components:
@@ -89,50 +98,65 @@
 #' @examples
 #'
 #' x0 <- c(1, 1)
-#' fn <- function(x) (x[1]-2)^2 + (x[2]-1)^2
-#' hin <- function(x) -0.25*x[1]^2 - x[2]^2 + 1  # hin >= 0
-#' heq <- function(x) x[1] - 2*x[2] + 1      # heq == 0
+#' fn <- function(x) (x[1] - 2) ^ 2 + (x[2] - 1) ^ 2
+#' hin <- function(x) 0.25 * x[1]^2 + x[2] ^ 2 - 1  # hin <= 0
+#' heq <- function(x) x[1] - 2 * x[2] + 1           # heq = 0
 #' gr <- function(x) nl.grad(x, fn)
 #' hinjac <- function(x) nl.jacobian(x, hin)
 #' heqjac <- function(x) nl.jacobian(x, heq)
 #'
-#' auglag(x0, fn, gr = NULL, hin = hin, heq = heq) # with COBYLA
+#' # with COBYLA
+#' auglag(x0, fn, gr = NULL, hin = hin, heq = heq, deprecatedBehavior = FALSE)
+#'
 #' # $par:   0.8228761 0.9114382
 #' # $value:   1.393464
 #' # $iter:  1001
 #'
-#' auglag(x0, fn, gr = NULL, hin = hin, heq = heq, localsolver = "SLSQP")
+#' auglag(x0, fn, gr = NULL, hin = hin, heq = heq, localsolver = "SLSQP",
+#'        deprecatedBehavior = FALSE)
+#'
 #' # $par:   0.8228757 0.9114378
 #' # $value:   1.393465
-#' # $iter   173
+#' # $iter   184
 #'
 #' ##  Example from the alabama::auglag help page
-#' fn <- function(x) (x[1] + 3*x[2] + x[3])^2 + 4 * (x[1] - x[2])^2
+#' ##  Parameters should be roughly (0, 0, 1) with an objective value of 1.
+#'
+#' fn <- function(x) (x[1] + 3 * x[2] + x[3]) ^ 2 + 4 * (x[1] - x[2]) ^ 2
 #' heq <- function(x) x[1] + x[2] + x[3] - 1
-#' hin <- function(x) c(6*x[2] + 4*x[3] - x[1]^3 - 3, x[1], x[2], x[3])
+#' # hin restated from alabama example to be <= 0.
+#' hin <- function(x) c(-6 * x[2] - 4 * x[3] + x[1] ^ 3 + 3, -x[1], -x[2], -x[3])
+#'
+#' set.seed(12)
+#' auglag(runif(3), fn, hin = hin, heq = heq, localsolver= "lbfgs",
+#'        deprecatedBehavior = FALSE)
 #'
-#' auglag(runif(3), fn, hin = hin, heq = heq, localsolver="lbfgs")
-#' # $par:   2.380000e-09 1.086082e-14 1.000000e+00
+#' # $par:   4.861756e-08 4.732373e-08 9.999999e-01
 #' # $value:   1
-#' # $iter:  289
+#' # $iter:  145
 #'
 #' ##  Powell problem from the Rsolnp::solnp help page
+#' ##  Parameters should be roughly (-1.7171, 1.5957, 1.8272, -0.7636, -0.7636)
+#' ##  with an objective value of 0.0539498478.
+#'
 #' x0 <- c(-2, 2, 2, -1, -1)
-#' fn1  <- function(x) exp(x[1]*x[2]*x[3]*x[4]*x[5])
+#' fn1  <- function(x) exp(x[1] * x[2] * x[3] * x[4] * x[5])
 #' eqn1 <-function(x)
-#' 	c(x[1]*x[1]+x[2]*x[2]+x[3]*x[3]+x[4]*x[4]+x[5]*x[5],
-#' 	  x[2]*x[3]-5*x[4]*x[5],
-#' 	  x[1]*x[1]*x[1]+x[2]*x[2]*x[2])
+#' 	c(x[1] * x[1] + x[2] * x[2] + x[3] * x[3] + x[4] * x[4] + x[5] * x[5] - 10,
+#' 	  x[2] * x[3] - 5 * x[4] * x[5],
+#' 	  x[1] * x[1] * x[1] + x[2] * x[2] * x[2] + 1)
 #'
-#' auglag(x0, fn1, heq = eqn1, localsolver = "mma")
-#' # $par: -3.988458e-10 -1.654201e-08 -3.752028e-10  8.904445e-10  8.926336e-10
-#' # $value:   1
-#' # $iter:  1001
+#' auglag(x0, fn1, heq = eqn1, localsolver = "mma", deprecatedBehavior = FALSE)
+#'
+#' # $par: -1.7173645  1.5959655  1.8268352 -0.7636185 -0.7636185
+#' # $value:   0.05394987
+#' # $iter:  916
 #'
 auglag <- function(x0, fn, gr = NULL, lower = NULL, upper = NULL, hin = NULL,
                    hinjac = NULL, heq = NULL, heqjac = NULL,
                    localsolver = "COBYLA", localtol = 1e-6, ineq2local = FALSE,
-                   nl.info = FALSE, control = list(), ...) {
+                   nl.info = FALSE, control = list(), deprecatedBehavior = TRUE,
+                   ...) {
   if (ineq2local) {
     # gsolver <- "NLOPT_LN_AUGLAG_EQ"
     stop("Inequalities to local solver: feature not yet implemented.")
@@ -167,23 +191,25 @@ auglag <- function(x0, fn, gr = NULL, lower = NULL, upper = NULL, hin = NULL,
 
   # Inequality constraints
   if (!is.null(hin)) {
-    if (getOption("nloptr.show.inequality.warning")) {
-      message("For consistency with the rest of the package the ",
-              "inequality sign may be switched from >= to <= in a ",
-              "future nloptr version.")
+    if (deprecatedBehavior) {
+      warning("The old behavior for hin >= 0 has been deprecated. Please ",
+              "restate the inequality to be <=0. The ability to use the old ",
+              "behavior will be removed in a future release.")
+      .hin <- match.fun(hin)
+      hin <- function(x) -.hin(x)      # change  hin >= 0  to  hin <= 0 !
     }
-
-    .hin <- match.fun(hin)
-    hin <- function(x) -.hin(x)   # change  hin >= 0  to  hin <= 0 !
   }
   if (!dfree) {
     if (is.null(hinjac)) {
       hinjac <- function(x) nl.jacobian(x, hin)
-    } else {
+    } else if (deprecatedBehavior) {
+      warning("The old behavior for hinjac >= 0 has been deprecated. Please ",
+              "restate the inequality to be <=0. The ability to use the old ",
+              "behavior will be removed in a future release.")
       .hinjac <- match.fun(hinjac)
       hinjac <- function(x) -.hinjac(x)
+      }
     }
-  }
 
   # Equality constraints
   if (!is.null(heq)) {

diff --git a/R/ccsaq.R b/R/ccsaq.R
@@ -6,7 +6,10 @@
 # Wrapper to solve optimization problem using CCSAQ.
 #
 # CHANGELOG:
-#   2023-02-11: Tweaks for efficiency and readability (Avraham Adler)
+# 2023-02-11: Tweaks for efficiency and readability (Avraham Adler)
+# 2024-06-04: Switched desired direction of the hin/hinjac inequalities, leaving
+#       the old behavior as the default for now. Also cleaned up the HS100
+#       example (Avraham Adler).
 #
 
 #' Conservative Convex Separable Approximation with Affine Approximation plus
@@ -28,6 +31,10 @@
 #' numerically if not specified.
 #' @param nl.info logical; shall the original NLopt info been shown.
 #' @param control list of options, see \code{nl.opts} for help.
+#' @param deprecatedBehavior logical; if \code{TRUE} (default for now), the old
+#' behavior of the Jacobian function is used, where the equality is \eqn{\ge 0}
+#' instead of \eqn{\le 0}. This will be reversed in a future release and
+#' eventually removed.
 #' @param ... additional arguments passed to the function.
 #'
 #' @return List with components:
@@ -54,48 +61,58 @@
 #' @examples
 #'
 #' ##  Solve the Hock-Schittkowski problem no. 100 with analytic gradients
+#' ##  See https://apmonitor.com/wiki/uploads/Apps/hs100.apm
+#'
 #' x0.hs100 <- c(1, 2, 0, 4, 0, 1, 1)
-#' fn.hs100 <- function(x) {
-#'   (x[1] - 10) ^ 2 + 5 * (x[2] - 12) ^ 2 + x[3] ^ 4 + 3 * (x[4] - 11) ^ 2 +
-#'    10 * x[5] ^ 6 + 7 * x[6] ^ 2 + x[7] ^ 4 - 4 * x[6] * x[7] -
-#'    10 * x[6] - 8 * x[7]
-#' }
-#' hin.hs100 <- function(x) {
-#'   h <- numeric(4)
-#'   h[1] <- 127 - 2 * x[1] ^ 2 - 3 * x[2] ^ 4 - x[3] - 4 * x[4] ^ 2 - 5 * x[5]
-#'   h[2] <- 282 - 7 * x[1] - 3 * x[2] - 10 * x[3] ^ 2 - x[4] + x[5]
-#'   h[3] <- 196 - 23 * x[1] - x[2] ^ 2 - 6 * x[6] ^ 2 + 8 * x[7]
-#'   h[4] <- -4 * x[1] ^ 2 - x[2] ^ 2 + 3 * x[1] * x[2] -2 * x[3] ^ 2 -
-#'        5 * x[6] + 11 * x[7]
-#'   return(h)
+#' fn.hs100 <- function(x) {(x[1] - 10) ^ 2 + 5 * (x[2] - 12) ^ 2 + x[3] ^ 4 +
+#'                          3 * (x[4] - 11) ^ 2 + 10 * x[5] ^ 6 + 7 * x[6] ^ 2 +
+#'                          x[7] ^ 4 - 4 * x[6] * x[7] - 10 * x[6] - 8 * x[7]}
+#'
+#' hin.hs100 <- function(x) {c(
+#' 2 * x[1] ^ 2 + 3 * x[2] ^ 4 + x[3] + 4 * x[4] ^ 2 + 5 * x[5] - 127,
+#' 7 * x[1] + 3 * x[2] + 10 * x[3] ^ 2 + x[4] - x[5] - 282,
+#' 23 * x[1] + x[2] ^ 2 + 6 * x[6] ^ 2 - 8 * x[7] - 196,
+#' 4 * x[1] ^ 2 + x[2] ^ 2 - 3 * x[1] * x[2] + 2 * x[3] ^ 2 + 5 * x[6] -
+#'  11 * x[7])
 #' }
+#'
 #' gr.hs100 <- function(x) {
-#'  c( 2 * x[1] -  20,
+#'  c( 2 * x[1] - 20,
 #'    10 * x[2] - 120,
 #'     4 * x[3] ^ 3,
 #'     6 * x[4] - 66,
 #'    60 * x[5] ^ 5,
 #'    14 * x[6] - 4 * x[7] - 10,
-#'     4 * x[7] ^ 3 - 4 * x[6] -  8)
+#'     4 * x[7] ^ 3 - 4 * x[6] - 8)
 #' }
+#'
 #' hinjac.hs100 <- function(x) {
-#'   matrix(c(4 * x[1], 12 * x[2] ^ 3, 1, 8 * x[4], 5, 0, 0, 7, 3, 20 * x[3],
-#'        1, -1, 0, 0, 23, 2 * x[2], 0, 0, 0, 12 * x[6], -8,
-#'        8 * x[1] - 3 * x[2], 2 * x[2] - 3 * x[1], 4 * x[3],
-#'        0, 0, 5, -11), 4, 7, byrow = TRUE)
+#'   matrix(c(4 * x[1], 12 * x[2] ^ 3, 1, 8 * x[4], 5, 0, 0,
+#'            7, 3, 20 * x[3], 1, -1, 0, 0,
+#'            23, 2 * x[2], 0, 0, 0, 12 * x[6], -8,
+#'            8 * x[1] - 3 * x[2], 2 * x[2] - 3 * x[1], 4 * x[3], 0, 0, 5, -11),
+#'            nrow = 4, byrow = TRUE)
 #' }
 #'
-#' # incorrect result with exact jacobian
+#' ##  The optimum value of the objective function should be 680.6300573
+#' ##  A suitable parameter vector is roughly
+#' ##  (2.330, 1.9514, -0.4775, 4.3657, -0.6245, 1.0381, 1.5942)
+#'
+#' # Results with exact Jacobian
 #' S <- ccsaq(x0.hs100, fn.hs100, gr = gr.hs100,
 #'       hin = hin.hs100, hinjac = hinjac.hs100,
-#'       nl.info = TRUE, control = list(xtol_rel = 1e-8))
+#'       nl.info = TRUE, control = list(xtol_rel = 1e-8),
+#'       deprecatedBehavior = FALSE)
 #'
-#' \donttest{
+#' # Results without Jacobian
 #' S <- ccsaq(x0.hs100, fn.hs100, hin = hin.hs100,
-#'       nl.info = TRUE, control = list(xtol_rel = 1e-8))
-#' }
+#'       nl.info = TRUE, control = list(xtol_rel = 1e-8),
+#'       deprecatedBehavior = FALSE)
+#'
+
 ccsaq <- function(x0, fn, gr = NULL, lower = NULL, upper = NULL, hin = NULL,
-                  hinjac = NULL, nl.info = FALSE, control = list(), ...) {
+                  hinjac = NULL, nl.info = FALSE, control = list(),
+                  deprecatedBehavior = TRUE, ...) {
 
   opts <- nl.opts(control)
   opts["algorithm"] <- "NLOPT_LD_CCSAQ"
@@ -111,19 +128,21 @@ ccsaq <- function(x0, fn, gr = NULL, lower = NULL, upper = NULL, hin = NULL,
   }
 
   if (!is.null(hin)) {
-    if (getOption("nloptr.show.inequality.warning")) {
-      message("For consistency with the rest of the package the ",
-              "inequality sign may be switched from >= to <= in a ",
-              "future nloptr version.")
+    if (deprecatedBehavior) {
+      warning("The old behavior for hin >= 0 has been deprecated. Please ",
+              "restate the inequality to be <=0. The ability to use the old ",
+              "behavior will be removed in a future release.")
+      .hin <- match.fun(hin)
+      hin <- function(x) -.hin(x)      # change  hin >= 0  to  hin <= 0 !
     }
-
-    .hin <- match.fun(hin)
-    hin <- function(x) -.hin(x)       # change  hin >= 0  to  hin <= 0 !
     if (is.null(hinjac)) {
       hinjac <- function(x) nl.jacobian(x, hin)
-    } else {
-      .hinjac <- match.fun(hinjac)
-      hinjac <- function(x) -.hinjac(x)
+    } else if (deprecatedBehavior) {
+      warning("The old behavior for hinjac >= 0 has been deprecated. Please ",
+              "restate the inequality to be <=0. The ability to use the old ",
+              "behavior will be removed in a future release.")
+        .hinjac <- match.fun(hinjac)
+        hinjac <- function(x) -.hinjac(x)
     }
   }