#'
#' The Constraint class.
#'
#' This virtual class represents a mathematical constraint.
#'
#' @name Constraint-class
#' @aliases Constraint
#' @rdname Constraint-class
setClass("Constraint", representation(dual_variables = "list"), prototype(dual_variables = list()), contains = "Canonical")
setMethod("initialize", "Constraint", function(.Object, ..., dual_variables = list()) {
.Object <- callNextMethod(.Object, ...)
# .Object@dual_variables <- lapply(.Object@args, function(arg) { Variable(dim(arg)) })
.Object@dual_variables <- lapply(.Object@args, function(arg) { new("Variable", dim = dim(arg)) })
return(.Object)
})
#' @param x,object A \linkS4class{Constraint} object.
#' @rdname Constraint-class
setMethod("as.character", "Constraint", function(x) { name(x) })
setMethod("show", "Constraint", function(object) {
print(paste(class(object), "(", as.character(object@args[[1]]), ")", sep = ""))
})
#' @describeIn Constraint The dimensions of the constrained expression.
setMethod("dim", "Constraint", function(x) { dim(x@args[[1]]) })
#' @describeIn Constraint The size of the constrained expression.
setMethod("size", "Constraint", function(object) { size(object@args[[1]]) })
#' @describeIn Constraint Is the constraint real?
setMethod("is_real", "Constraint", function(object) { !is_complex(object) })
#' @describeIn Constraint Is the constraint imaginary?
setMethod("is_imag", "Constraint", function(object) { all(sapply(object@args, is_imag)) })
#' @describeIn Constraint Is the constraint complex?
setMethod("is_complex", "Constraint", function(object) { any(sapply(object@args, is_complex)) })
#' @describeIn Constraint Is the constraint DCP?
setMethod("is_dcp", "Constraint", function(object) { stop("Unimplemented") })
#' @describeIn Constraint Is the constraint DGP?
setMethod("is_dgp", "Constraint", function(object) { stop("Unimplemented") })
#' @describeIn Constraint The residual of a constraint
setMethod("residual", "Constraint", function(object) { stop("Unimplemented") })
#' @describeIn Constraint The violation of a constraint.
setMethod("violation", "Constraint", function(object) {
resid <- residual(object)
if(any(is.na(resid)))
stop("Cannot compute the violation of a constraint whose expression is NA-valued.")
return(resid)
})
#' @param tolerance The tolerance for checking if the constraint is violated.
#' @describeIn Constraint The value of a constraint.
setMethod("constr_value", "Constraint", function(object, tolerance = 1e-8) {
resid <- residual(object)
if(any(is.na(resid)))
stop("Cannot compute the value of a constraint whose expression is NA-valued.")
return(all(resid <= tolerance))
})
#'
#' A Class Union of List and Constraint
#'
#' @name ListORConstr-class
#' @rdname ListORConstr-class
setClassUnion("ListORConstr", c("list", "Constraint"))
# Helper function since syntax is different for LinOp (list) vs. Constraint object
#' @param object A list or \linkS4class{Constraint} object.
#' @describeIn ListORConstr Returns the ID associated with the list or constraint.
setMethod("id", "ListORConstr", function(object) {
if(is.list(object))
object$constr_id
else
object@id
})
#' @describeIn Constraint Information needed to reconstruct the object aside from the args.
setMethod("get_data", "Constraint", function(object) { list(id(object)) })
#' @describeIn Constraint The dual values of a constraint.
setMethod("dual_value", "Constraint", function(object) { value(object@dual_variables[[1]]) })
#' @param value A numeric scalar, vector, or matrix.
#' @describeIn Constraint Replaces the dual values of a constraint..
setReplaceMethod("dual_value", "Constraint", function(object, value) {
object@dual_variables[[1]] <- value
object
})
#'
#' The ZeroConstraint class
#'
#' @rdname ZeroConstraint-class
.ZeroConstraint <- setClass("ZeroConstraint", representation(expr = "Expression"), contains = "Constraint")
ZeroConstraint <- function(expr, id = NA_integer_) { .ZeroConstraint(expr = expr, id = id) }
setMethod("initialize", "ZeroConstraint", function(.Object, ..., expr) {
.Object@expr <- expr
callNextMethod(.Object, ..., args = list(expr))
})
#' @param x,object A \linkS4class{ZeroConstraint} object.
#' @describeIn ZeroConstraint The string representation of the constraint.
setMethod("name", "ZeroConstraint", function(x) {
# paste(as.character(x@args[[1]]), "== 0")
paste(name(x@args[[1]]), "== 0")
})
#' @describeIn ZeroConstraint The dimensions of the constrained expression.
setMethod("dim", "ZeroConstraint", function(x) { dim(x@args[[1]]) })
#' @describeIn Constraint The size of the constrained expression.
setMethod("size", "ZeroConstraint", function(object) { size(object@args[[1]]) })
#' @describeIn ZeroConstraint Is the constraint DCP?
setMethod("is_dcp", "ZeroConstraint", function(object) { is_affine(object@args[[1]]) })
#' @describeIn ZeroConstraint Is the constraint DGP?
setMethod("is_dgp", "ZeroConstraint", function(object) { FALSE })
#' @describeIn ZeroConstraint The residual of a constraint
setMethod("residual", "ZeroConstraint", function(object) {
val <- value(expr(object))
if(any(is.na(val)))
return(NA_real_)
return(abs(val))
})
#' @describeIn ZeroConstraint The graph implementation of the object.
setMethod("canonicalize", "ZeroConstraint", function(object) {
canon <- canonical_form(object@args[[1]])
obj <- canon[[1]]
constraints <- canon[[2]]
dual_holder <- create_eq(obj, constr_id = id(object))
return(list(NA, c(constraints, list(dual_holder))))
})
#'
#' The EqConstraint class
#'
#' @rdname EqConstraint-class
.EqConstraint <- setClass("EqConstraint", representation(lhs = "ConstValORExpr", rhs = "ConstValORExpr", expr = "ConstValORExpr"), prototype(expr = NA_real_), contains = "Constraint")
EqConstraint <- function(lhs, rhs, id = NA_integer_) { .EqConstraint(lhs = lhs, rhs = rhs, id = id) }
setMethod("initialize", "EqConstraint", function(.Object, ..., lhs, rhs, expr = NA_real_) {
.Object@lhs <- lhs
.Object@rhs <- rhs
.Object@expr <- lhs - rhs
callNextMethod(.Object, ..., args = list(lhs, rhs))
})
setMethod(".construct_dual_variables", "EqConstraint", function(object, args) {
callNextMethod(object, list(object@expr))
})
#' @param x,object A \linkS4class{EqConstraint} object.
#' @describeIn EqConstraint The string representation of the constraint.
setMethod("name", "EqConstraint", function(x) {
# paste(as.character(x@args[[1]]), "==", as.character(x@args[[2]]))
paste(name(x@args[[1]]), "==", name(x@args[[2]]))
})
#' @describeIn EqConstraint The dimensions of the constrained expression.
setMethod("dim", "EqConstraint", function(x) { dim(x@expr) })
#' @describeIn EqConstraint The size of the constrained expression.
setMethod("size", "EqConstraint", function(object) { size(object@expr) })
#' @describeIn EqConstraint The expression to constrain.
setMethod("expr", "EqConstraint", function(object) { object@expr })
#' @describeIn EqConstraint Is the constraint DCP?
setMethod("is_dcp", "EqConstraint", function(object) { is_affine(object@expr) })
#' @describeIn EqConstraint Is the constraint DGP?
setMethod("is_dgp", "EqConstraint", function(object) {
is_log_log_affine(object@args[[1]]) && is_log_log_affine(object@args[[2]])
})
#' @describeIn EqConstraint The residual of the constraint..
setMethod("residual", "EqConstraint", function(object) {
val <- value(object@expr)
if(any(is.na(val)))
return(NA_real_)
return(abs(val))
})
#'
#' The NonPosConstraint class
#'
#' @rdname NonPosConstraint-class
.NonPosConstraint <- setClass("NonPosConstraint", representation(expr = "Expression"), contains = "Constraint")
NonPosConstraint <- function(expr, id = NA_integer_) { .NonPosConstraint(expr = expr, id = id) }
setMethod("initialize", "NonPosConstraint", function(.Object, ..., expr) {
.Object@expr <- expr
callNextMethod(.Object, ..., args = list(expr))
})
#' @param x,object A \linkS4class{NonPosConstraint} object.
#' @describeIn NonPosConstraint The string representation of the constraint.
setMethod("name", "NonPosConstraint", function(x) {
# paste(as.character(x@args[[1]]), "<= 0")
paste(name(x@args[[1]]), "<= 0")
})
#' @describeIn NonPosConstraint Is the constraint DCP?
setMethod("is_dcp", "NonPosConstraint", function(object) { is_convex(object@args[[1]]) })
#' @describeIn NonPosConstraint Is the constraint DGP?
setMethod("is_dgp", "NonPosConstraint", function(object) { FALSE })
#' @describeIn NonPosConstraint The graph implementation of the object.
setMethod("canonicalize", "NonPosConstraint", function(object) {
canon <- canonical_form(object@args[[1]])
obj <- canon[[1]]
constraints <- canon[[2]]
dual_holder <- create_leq(obj, constr_id = id(object))
return(list(NA, c(constraints, list(dual_holder))))
})
#' @describeIn NonPosConstraint The residual of the constraint.
setMethod("residual", "NonPosConstraint", function(object) {
val <- value(expr(object))
if(any(is.na(val)))
return(NA_real_)
return(pmax(val, 0))
})
#'
#' The IneqConstraint class
#'
#' @rdname IneqConstraint-class
.IneqConstraint <- setClass("IneqConstraint", representation(lhs = "ConstValORExpr", rhs = "ConstValORExpr", expr = "ConstValORExpr"), prototype(expr = NA_real_), contains = "Constraint")
IneqConstraint <- function(lhs, rhs, id = NA_integer_) { .IneqConstraint(lhs = lhs, rhs = rhs, id = id) }
setMethod("initialize", "IneqConstraint", function(.Object, ..., lhs, rhs, expr = NA_real_) {
.Object@lhs <- lhs
.Object@rhs <- rhs
.Object@expr <- lhs - rhs
if(is_complex(.Object@expr))
stop("Inequality constraints cannot be complex.")
callNextMethod(.Object, ..., args = list(lhs, rhs))
})
#' @param x,object A \linkS4class{IneqConstraint} object.
#' @describeIn IneqConstraint The string representation of the constraint.
setMethod("name", "IneqConstraint", function(x) {
# paste(as.character(x@args[[1]]), "<=", as.character(x@args[[2]]))
paste(name(x@args[[1]]), "<=", name(x@args[[2]]))
})
#' @describeIn IneqConstraint The dimensions of the constrained expression.
setMethod("dim", "IneqConstraint", function(x) { dim(x@expr) })
#' @describeIn IneqConstraint The size of the constrained expression.
setMethod("size", "IneqConstraint", function(object) { size(object@expr) })
#' @describeIn IneqConstraint The expression to constrain.
setMethod("expr", "IneqConstraint", function(object) { object@expr })
#' @describeIn IneqConstraint A non-positive constraint is DCP if its argument is convex.
setMethod("is_dcp", "IneqConstraint", function(object) { is_convex(object@expr) })
#' @describeIn IneqConstraint Is the constraint DGP?
setMethod("is_dgp", "IneqConstraint", function(object) {
is_log_log_convex(object@args[[1]]) && is_log_log_concave(object@args[[2]])
})
#' @describeIn IneqConstraint The residual of the constraint.
setMethod("residual", "IneqConstraint", function(object) {
val <- value(object@expr)
if(any(is.na(val)))
return(NA_real_)
return(pmax(val, 0))
})
# TODO: Do I need the NonlinearConstraint class?
#'
#' The NonlinearConstraint class.
#'
#' This class represents a nonlinear inequality constraint, \eqn{f(x) \leq 0} where \eqn{f} is twice-differentiable.
#'
#' @slot f A nonlinear function.
#' @slot vars_ A list of variables involved in the function.
#' @slot .x_dim (Internal) The dimensions of a column vector with number of elements equal to the total elements in all the variables.
#' @name NonlinearConstraint-class
#' @aliases NonlinearConstraint
#' @rdname NonlinearConstraint-class
.NonlinearConstraint <- setClass("NonlinearConstraint", representation(f = "function", vars_ = "list", .x_dim = "numeric"),
prototype(.x_dim = NULL), contains = "Constraint")
#' @param f A nonlinear function.
#' @param vars_ A list of variables involved in the function.
#' @param id (Optional) An integer representing the unique ID of the contraint.
#' @rdname NonlinearConstraint-class
NonlinearConstraint <- function(f, vars_, id = NA_integer_) { .NonlinearConstraint(f = f, vars_ = vars_, id = id) }
setMethod("initialize", "NonlinearConstraint", function(.Object, ..., f, vars_) {
.Object@f <- f
.Object@vars_ <- vars_
# The dimensions of vars_ in f(vars_)
sizes <- sapply(.Object@vars_, function(v) { as.integer(prod(dim(v))) })
.Object@.x_dim <- c(sum(sizes), 1)
callNextMethod(.Object, ..., args = .Object@vars_)
})
# Add the block to a slice of the matrix.
setMethod("block_add", "NonlinearConstraint", function(object, mat, block, vert_offset, horiz_offset, rows, cols, vert_step = 1, horiz_step = 1) {
if(is(mat, "sparseMatrix") && is.matrix(block))
block <- Matrix(block, sparse = TRUE)
row_seq <- seq(vert_offset, rows + vert_offset, vert_step)
col_seq <- seq(horiz_offset, cols + horiz_offset, horiz_step)
mat[row_seq, col_seq] <- mat[row_seq, col_seq] + block
mat
})
# Place \code{x_0 = f()} in the vector of all variables.
setMethod("place_x0", "NonlinearConstraint", function(object, big_x, var_offsets) {
tmp <- object@f()
m <- tmp[[1]]
x0 <- tmp[[2]]
offset <- 0
for(var in object@args) {
var_dim <- as.integer(prod(dim(var)))
var_x0 <- x0[offset:(offset + var_dim)]
big_x <- block_add(object, big_x, var_x0, var_offsets[get_data(var)], 0, var_dim, 1)
offset <- offset + var_dim
}
big_x
})
# Place \code{Df} in the gradient of all functions.
setMethod("place_Df", "NonlinearConstraint", function(object, big_Df, Df, var_offsets, vert_offset) {
horiz_offset <- 0
for(var in object@args) {
var_dim <- as.integer(prod(dim(var)))
var_Df <- Df[, horiz_offset:(horiz_offset + var_dim)]
big_Df <- block_add(object, big_Df, var_Df, vert_offset, var_offsets[get_data(var)], num_cones(object), var_dim)
horiz_offset <- horiz_offset + var_dim
}
big_Df
})
# Place \code{H} in the Hessian of all functions.
setMethod("place_H", "NonlinearConstraint", function(object, big_H, H, var_offsets) {
offset <- 0
for(var in object@args) {
var_dim <- as.integer(prod(dim(var)))
var_H <- H[offset:(offset + var_dim), offset:(offset + var_dim)]
big_H <- block_add(object, big_H, var_H, var_offsets[get_data(var)], var_offsets[get_data(var)], var_dim, var_dim)
offset <- offset + var_dim
}
big_H
})
# Extract the function variables from the vector \code{x} of all variables.
setMethod("extract_variables", "NonlinearConstraint", function(object, x, var_offsets) {
local_x <- matrix(0, nrow = object@.x_dim[1], ncol = object@.x_dim[2])
offset <- 0
for(var in object@args) {
var_dim <- as.integer(prod(dim(var)))
value <- x[var_offsets[get_data(var)]:(var_offsets[get_data(var)] + var_dim)]
local_x <- block_add(object, local_x, value, offset, 0, var_dim, 1)
offset <- offset + var_dim
}
local_x
})
#'
#' The ExpCone class.
#'
#' This class represents a reformulated exponential cone constraint operating elementwise on \eqn{a, b, c}.
#'
#' Original cone:
#' \deqn{
#' K = \{(x,y,z) | y > 0, ye^{x/y} \leq z\} \cup \{(x,y,z) | x \leq 0, y = 0, z \geq 0\}
#' }
#' Reformulated cone:
#' \deqn{
#' K = \{(x,y,z) | y, z > 0, y\log(y) + x \leq y\log(z)\} \cup \{(x,y,z) | x \leq 0, y = 0, z \geq 0\}
#' }
#'
#' @slot x The variable \eqn{x} in the exponential cone.
#' @slot y The variable \eqn{y} in the exponential cone.
#' @slot z The variable \eqn{z} in the exponential cone.
#' @name ExpCone-class
#' @aliases ExpCone
#' @rdname ExpCone-class
.ExpCone <- setClass("ExpCone", representation(x = "ConstValORExpr", y = "ConstValORExpr", z = "ConstValORExpr"), contains = "Constraint")
#' @param x The variable \eqn{x} in the exponential cone.
#' @param y The variable \eqn{y} in the exponential cone.
#' @param z The variable \eqn{z} in the exponential cone.
#' @param id (Optional) A numeric value representing the constraint ID.
#' @rdname ExpCone-class
## #' @export
ExpCone <- function(x, y, z, id = NA_integer_) { .ExpCone(x = x, y = y, z = z, id = id) }
setMethod("initialize", "ExpCone", function(.Object, ..., x, y, z) {
.Object@x <- x
.Object@y <- y
.Object@z <- z
callNextMethod(.Object, ..., args = list(.Object@x, .Object@y, .Object@z))
})
setMethod("show", "ExpCone", function(object) {
print(paste("ExpCone(", as.character(object@x), ", ", as.character(object@y), ", ", as.character(object@z), ")", sep = ""))
})
#' @rdname ExpCone-class
setMethod("as.character", "ExpCone", function(x) {
paste("ExpCone(", as.character(x@x), ", ", as.character(x@y), ", ", as.character(x@z), ")", sep = "")
})
#' @param object A \linkS4class{ExpCone} object.
#' @describeIn ExpCone The size of the \code{x} argument.
setMethod("residual", "ExpCone", function(object) {
# TODO: The projection should be implemented directly.
if(any(is.na(value(object@x))) || any(is.na(value(object@y))) || any(is.na(value(object@z))))
return(NA_real_)
# x <- Variable(dim(object@x))
# y <- Variable(dim(object@y))
# z <- Variable(dim(object@z))
x <- new("Variable", dim = dim(object@x))
y <- new("Variable", dim = dim(object@y))
z <- new("Variable", dim = dim(object@z))
constr <- list(ExpCone(x, y, z))
obj <- Minimize(Norm2(HStack(x, y, z) - HStack(value(object@x), value(object@y), value(object@z))))
prob <- Problem(obj, constr)
result <- solve(prob)
return(result$value)
})
#' @describeIn ExpCone The number of entries in the combined cones.
setMethod("size", "ExpCone", function(object) {
# TODO: Use size of dual variable(s) instead.
sum(cone_sizes(object))
})
#' @describeIn ExpCone The number of elementwise cones.
setMethod("num_cones", "ExpCone", function(object) {
as.integer(prod(dim(object@args[[1]])))
})
#' @describeIn ExpCone The dimensions of the exponential cones.
setMethod("cone_sizes", "ExpCone", function(object) {
rep(num_cones(object), 3)
})
#' @describeIn ExpCone An exponential constraint is DCP if each argument is affine.
setMethod("is_dcp", "ExpCone", function(object) {
all(sapply(object@args, is_affine))
})
#' @describeIn ExpCone Is the constraint DGP?
setMethod("is_dgp", "ExpCone", function(object) { FALSE })
#' @describeIn ExpCone Canonicalizes by converting expressions to LinOps.
setMethod("canonicalize", "ExpCone", function(object) {
arg_objs <- list()
arg_constr <- list()
for(arg in object@args) {
canon <- canonical_form(arg)
arg_objs <- c(arg_objs, canon[[1]])
arg_constr <- c(arg_constr, canon[[2]])
}
exp_constr <- do.call(ExpCone, arg_objs)
list(0, c(list(exp_constr), arg_constr))
})
#'
#' The PSDConstraint class.
#'
#' This class represents the positive semidefinite constraint, \eqn{\frac{1}{2}(X + X^T) \succeq 0}, i.e. \eqn{z^T(X + X^T)z \geq 0} for all \eqn{z}.
#'
#' @slot expr An \linkS4class{Expression}, numeric element, vector, or matrix representing \eqn{X}.
#' @name PSDConstraint-class
#' @aliases PSDConstraint
#' @rdname PSDConstraint-class
.PSDConstraint <- setClass("PSDConstraint", representation(expr = "ConstValORExpr"),
validity = function(object) {
expr_dim <- dim(object@expr)
if(length(expr_dim) != 2 || expr_dim[1] != expr_dim[2])
stop("Non-square matrix in positive definite constraint.")
return(TRUE)
}, contains = "Constraint")
#' @param expr An \linkS4class{Expression}, numeric element, vector, or matrix representing \eqn{X}.
#' @param id (Optional) A numeric value representing the constraint ID.
#' @rdname PSDConstraint-class
PSDConstraint <- function(expr, id = NA_integer_) { .PSDConstraint(expr = expr, id = id) }
setMethod("initialize", "PSDConstraint", function(.Object, ..., expr) {
.Object@expr <- expr
callNextMethod(.Object, ..., args = list(expr))
})
#' @param x,object A \linkS4class{PSDConstraint} object.
#' @describeIn PSDConstraint The string representation of the constraint.
setMethod("name", "PSDConstraint", function(x) {
# paste(as.character(x@args[[1]]), ">> 0")
paste(name(x@args[[1]]), ">> 0")
})
#' @describeIn PSDConstraint The constraint is DCP if the left-hand and right-hand expressions are affine.
setMethod("is_dcp", "PSDConstraint", function(object) { is_affine(object@args[[1]]) })
#' @describeIn PSDConstraint Is the constraint DGP?
setMethod("is_dgp", "PSDConstraint", function(object) { FALSE })
#' @describeIn PSDConstraint A \linkS4class{Expression} representing the residual of the constraint.
setMethod("residual", "PSDConstraint", function(object) {
val <- value(expr(object))
if(any(is.na(val)))
return(NA_real_)
min_eig <- LambdaMin(object@args[[1]] + t(object@args[[1]]))/2
value(Neg(min_eig))
})
#' @describeIn PSDConstraint The graph implementation of the object. Marks the top level constraint as the \code{dual_holder} so the dual value will be saved to the \linkS4class{PSDConstraint}.
setMethod("canonicalize", "PSDConstraint", function(object) {
canon <- canonical_form(object@args[[1]])
obj <- canon[[1]]
constraints <- canon[[2]]
dual_holder <- PSDConstraint(obj, id = id(object))
return(list(NA, c(constraints, list(dual_holder))))
})
setMethod("format_constr", "PSDConstraint", function(object, eq_constr, leq_constr, dims, solver) {
.format <- function(object) {
leq_constr <- create_geq(expr(object), constr_id = object@id)
return(list(leq_constr))
}
new_leq_constr <- .format(object)
# 0 <= A.
leq_constr <- c(leq_constr, new_leq_constr)
# Update dims.
dims[[PSD_DIM]] <- c(dims[[PSD_DIM]], nrow(object))
list(eq_constr = eq_constr, leq_constr = leq_constr, dims = dims)
})
#'
#' The SOC class.
#'
#' This class represents a second-order cone constraint, i.e. \eqn{\|x\|_2 \leq t}.
#'
#' @slot t The scalar part of the second-order constraint.
#' @slot X A matrix whose rows/columns are each a cone.
#' @slot axis The dimension along which to slice: \code{1} indicates rows, and \code{2} indicates columns. The default is \code{2}.
#' @name SOC-class
#' @aliases SOC
#' @rdname SOC-class
.SOC <- setClass("SOC", representation(t = "ConstValORExpr", X = "ConstValORExpr", axis = "numeric"),
prototype(t = NA_real_, X = NA_real_, axis = 2), contains = "Constraint")
#' @param t The scalar part of the second-order constraint.
#' @param X A matrix whose rows/columns are each a cone.
#' @param axis The dimension along which to slice: \code{1} indicates rows, and \code{2} indicates columns. The default is \code{2}.
#' @param id (Optional) A numeric value representing the constraint ID.
#' @rdname SOC-class
## #' @export
SOC <- function(t, X, axis = 2, id = NA_integer_) { .SOC(t = t, X = X, axis = axis, id = id) }
setMethod("initialize", "SOC", function(.Object, ..., t, X, axis = 2) {
# TODO: Allow imaginary X.
# if(!(is.null(dim(t)) || length(dim(t)) == 1))
t_dim <- dim(t)
if(!(is.null(t_dim) || length(t_dim) == 1 || (length(t_dim) == 2 && t_dim[2] == 1)))
stop("t must be a scalar or 1-dimensional vector.")
.Object@t <- t
.Object@X <- X
.Object@axis <- axis
callNextMethod(.Object, ..., args = list(t, X))
})
#' @param x,object A \linkS4class{SOC} object.
#' @rdname SOC-class
setMethod("as.character", "SOC", function(x) {
paste("SOC(", as.character(x@t), ", ", as.character(x@X), ")", sep = "")
})
#' @describeIn SOC The residual of the second-order constraint.
setMethod("residual", "SOC", function(object) {
t <- value(object@args[[1]])
X <- value(object@args[[2]])
if(is.na(t) || is.na(X))
return(NA)
if(object@axis == 2)
X <- t(X)
norms <- apply(X, 1, function(row) { norm(row, "2") })
zero_indices <- which(X <= -t)[1]
averaged_indices <- which(X >= abs(t))[1]
X_proj <- as.matrix(X)
t_proj <- as.matrix(t)
X_proj[zero_indices] <- 0
t_proj[zero_indices] <- 0
avg_coeff <- 0.5*(1 + t/norms)
X_proj[averaged_indices] <- avg_coeff * X[averaged_indices]
t_proj[averaged_indices] <- avg_coeff * t[averaged_indices]
Xt_diff <- cbind(X, t) - cbind(X_proj, t_proj)
apply(Xt_diff, 1, function(col) { norm(col, "2") })
})
#' @describeIn SOC Information needed to reconstruct the object aside from the args.
setMethod("get_data", "SOC", function(object) { list(object@axis) })
#' @param eq_constr A list of the equality constraints in the canonical problem.
#' @param leq_constr A list of the inequality constraints in the canonical problem.
#' @param dims A list with the dimensions of the conic constraints.
#' @param solver A string representing the solver to be called.
#' @describeIn SOC Format SOC constraints as inequalities for the solver.
setMethod("format_constr", "SOC", function(object, eq_constr, leq_constr, dims, solver) {
.format <- function(object) {
list(list(), format_axis(object@args[[1]], object@args[[2]], object@axis))
}
leq_constr <- c(leq_constr, .format(object)[[2]])
dims[[SOC_DIM]] <- c(dims[[SOC_DIM]], cone_sizes(object))
list(eq_constr = eq_constr, leq_constr = leq_constr, dims = dims)
})
#' @describeIn SOC The number of elementwise cones.
setMethod("num_cones", "SOC", function(object) { prod(dim(object@args[[1]])) })
#' @describeIn SOC The number of entries in the combined cones.
setMethod("size", "SOC", function(object) {
# TODO: Use size of dual variable(s) instead.
sum(cone_sizes(object))
})
#' @describeIn SOC The dimensions of the second-order cones.
setMethod("cone_sizes", "SOC", function(object) {
if(object@axis == 2) # Collapse columns.
idx <- 1
else if(object@axis == 1) # Collapse rows.
idx <- 2
else
stop("Unimplemented")
cone_size <- 1 + dim(object@args[[2]])[idx]
sapply(1:num_cones(object), function(i) { cone_size })
})
#' @describeIn SOC An SOC constraint is DCP if each of its arguments is affine.
setMethod("is_dcp", "SOC", function(object) {
all(sapply(object@args, function(arg) { is_affine(arg) }))
})
#' @describeIn SOC Is the constraint DGP?
setMethod("is_dgp", "SOC", function(object) { FALSE })
# TODO: Hack.
#' @describeIn SOC The canonicalization of the constraint.
setMethod("canonicalize", "SOC", function(object) {
canon_t <- canonical_form(object@args[[1]])
t <- canon_t[[1]]
t_cons <- canon_t[[2]]
canon_X <- canonical_form(object@args[[2]])
X <- canon_X[[1]]
X_cons <- canon_X[[2]]
new_soc <- SOC(t, X, object@axis)
return(list(NA, c(list(new_soc), t_cons, X_cons)))
})
#'
#' The SOCAxis class.
#'
#' This class represents a second-order cone constraint for each row/column.
#' It Assumes \eqn{t} is a vector the same length as \eqn{X}'s rows (columns) for axis == 1 (2).
#'
#' @slot t The scalar part of the second-order constraint.
#' @slot x_elems A list containing \code{X}, a matrix whose rows/columns are each a cone.
#' @slot axis The dimension across which to take the slice: \code{1} indicates rows, and \code{2} indicates columns.
#' @name SOCAxis-class
#' @aliases SOCAxis
#' @rdname SOCAxis-class
.SOCAxis <- setClass("SOCAxis", representation(x_elems = "list"), prototype(x_elems = list()), contains = "SOC")
#' @param t The scalar part of the second-order constraint.
#' @param X A matrix whose rows/columns are each a cone.
#' @param axis The dimension across which to take the slice: \code{1} indicates rows, and \code{2} indicates columns.
#' @param id (Optional) A numeric value representing the constraint ID.
#' @rdname SOCAxis-class
## #' @export
SOCAxis <- function(t, X, axis, id = NA_integer_) { .SOCAxis(t = t, X = X, axis = axis, id = id) }
setMethod("initialize", "SOCAxis", function(.Object, ...) {
.Object <- callNextMethod(.Object, ...)
.Object@x_elems <- list(.Object@X)
.Object
})
#' @param x,object A \linkS4class{SOCAxis} object.
#' @rdname SOCAxis-class
setMethod("as.character", "SOCAxis", function(x) {
paste("SOCAxis(", as.character(x@t), ", ", as.character(x@X), ", <", paste(x@axis, collapse = ", "), ">)", sep = "")
})
#' @param eq_constr A list of the equality constraints in the canonical problem.
#' @param leq_constr A list of the inequality constraints in the canonical problem.
#' @param dims A list with the dimensions of the conic constraints.
#' @param solver A string representing the solver to be called.
#' @describeIn SOCAxis Format SOC constraints as inequalities for the solver.
setMethod("format_constr", "SOCAxis", function(object, eq_constr, leq_constr, dims, solver) {
.format <- function(object) {
list(list(), format_axis(object@t, object@x_elems[[1]], object@axis))
}
leq_constr <- c(leq_constr, .format(object)[[2]])
# Update dims
for(cone_size in size(object))
dims[[SOC_DIM]] <- c(dims[[SOC_DIM]], cone_size[1])
list(eq_constr = eq_constr, leq_constr = leq_constr, dims = dims)
})
#' @describeIn SOCAxis The number of elementwise cones.
setMethod("num_cones", "SOCAxis", function(object) { nrow(object@t) })
#' @describeIn SOCAxis The dimensions of a single cone.
setMethod("cone_sizes", "SOCAxis", function(object) {
if(object@axis == 1) # Return ncols if applying along each row
c(1 + ncol(object@x_elems[[1]]), 1)
else if(object@axis == 2) # Return nrows if applying along each column
c(1 + nrow(object@x_elems[[1]]), 1)
else
stop("Invalid axis ", object@axis)
})
#' @describeIn SOCAxis The dimensions of the (elementwise) second-order cones.
setMethod("size", "SOCAxis", function(object) {
cone_size <- cone_size(object)
lapply(1:num_cones(object), function(i) { cone_size })
})
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.