-
Notifications
You must be signed in to change notification settings - Fork 25
/
add_cuts_portfolio.R
129 lines (127 loc) · 4.36 KB
/
add_cuts_portfolio.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
#' @include Portfolio-proto.R
NULL
#' Add Bender's cuts portfolio
#'
#' Generate a portfolio of solutions for a conservation planning
#' [problem()] using Bender's cuts (discussed in Rodrigues
#' *et al.* 2000). This is recommended as a replacement for
#' [add_gap_portfolio()] when the *Gurobi* software is not
#' available.
#'
#' @param x [problem()] (i.e., [`ConservationProblem-class`]) object.
#'
#' @param number_solutions `integer` number of attempts to generate
#' different solutions. Defaults to 10.
#'
#' @details This strategy for generating a portfolio of solutions involves
#' solving the problem multiple times and adding additional constraints
#' to forbid previously obtained solutions. In general, this strategy is most
#' useful when problems take a long time to solve and benefit from
#' having multiple threads allocated for solving an individual problem.
#
#' @section Notes:
#' In early versions (< 4.0.1), this function was only compatible with
#' *Gurobi* (i.e., [add_gurobi_solver()]). To provide functionality with
#' exact algorithm solvers, this function now adds constraints to the
#' problem formulation to generate multiple solutions.
#'
#' @return Object (i.e., [`ConservationProblem-class`]) with the portfolio
#' added to it.
#'
#' @seealso
#' See [portfolios] for an overview of all functions for adding a portfolio.
#'
#' @family portfolios
#'
#' @references
#' Rodrigues AS, Cerdeira OJ, and Gaston KJ (2000) Flexibility,
#' efficiency, and accountability: adapting reserve selection algorithms to
#' more complex conservation problems. *Ecography*, 23: 565--574.
#'
#' @examples
#' # set seed for reproducibility
#' set.seed(500)
#'
#' # load data
#' data(sim_pu_raster, sim_features, sim_pu_zones_stack, sim_features_zones)
#'
#' # create minimal problem with cuts portfolio
#' p1 <- problem(sim_pu_raster, sim_features) %>%
#' add_min_set_objective() %>%
#' add_relative_targets(0.2) %>%
#' add_cuts_portfolio(10) %>%
#' add_default_solver(gap = 0.2, verbose = FALSE)
#'
#' \dontrun{
#' # solve problem and generate 10 solutions within 20% of optimality
#' s1 <- solve(p1)
#'
#' # plot solutions in portfolio
#' plot(stack(s1), axes = FALSE, box = FALSE)
#' }
#' # build multi-zone conservation problem with cuts portfolio
#' p2 <- problem(sim_pu_zones_stack, sim_features_zones) %>%
#' add_min_set_objective() %>%
#' add_relative_targets(matrix(runif(15, 0.1, 0.2), nrow = 5,
#' ncol = 3)) %>%
#' add_binary_decisions() %>%
#' add_cuts_portfolio(10) %>%
#' add_default_solver(gap = 0.2, verbose = FALSE)
#'
#' \dontrun{
#' # solve the problem
#' s2 <- solve(p2)
#'
#' # print solution
#' str(s2, max.level = 1)
#'
#' # plot solutions in portfolio
#' plot(stack(lapply(s2, category_layer)), main = "solution", axes = FALSE,
#' box = FALSE)
#' }
#' @name add_cuts_portfolio
NULL
#' @rdname add_cuts_portfolio
#' @export
add_cuts_portfolio <- function(x, number_solutions = 10L) {
# assert that arguments are valid
assertthat::assert_that(inherits(x, "ConservationProblem"),
assertthat::is.count(number_solutions),
isTRUE(all(is.finite(number_solutions))))
# add portfolio
x$add_portfolio(pproto(
"CutsPortfolio",
Portfolio,
name = "Cuts portfolio",
parameters = parameters(
integer_parameter("number_solutions", number_solutions,
lower_limit = 1L)),
run = function(self, x, solver) {
## extract number of desired solutions
n <- self$parameters$get("number_solutions")
## solve problem to verify that it is feasible
sol <- solver$solve(x)
## if solving the problem failed then return NULL
if (is.null(sol))
return(sol)
## generate additional solutions
sol <- list(sol)
for (i in seq_len(n - 1) + 1) {
### add cuts
rcpp_forbid_solution(x$ptr, sol[[i - 1]][[1]])
### solve solution
curr_sol <- solver$solve(x)
### if contains valid solution then
if (!is.null(curr_sol$x)) {
sol[[i]] <- curr_sol
} else {
if ((i + 1) < self$parameters$get("number_solutions"))
warning(paste("there are only", length(sol),
"feasible solutions within the optimality gap"))
break()
}
}
## compile results
return(sol)
}
))}