D1 <- D2 <- diag(0, p)
D1[1:5, 1:5] <- 1; D2[5:p, 5:p] <- 1;
D3 <- D4 <- diag(0, p - 1)
D3[1:4, 1:4] <- 1; D4[5:(p-1), 5:(p-1)] <- 1;
constraints = list(D3 %*% diff(D1 %*% beta) >= 0, D4 %*% diff(D2 %*% beta) <= 0)
problem4 <- Problem(objective, constraints)
ex4 <- solve(problem4)
betaEstimate <- ex4$getValue(beta)

# generate (p-1) x p matrix for both constraints
# beta_1 <= beta_2 implies beta_2 - beta_1 >= 0, etc.
# and beta_5 >= beta_6 implies beta_5 - beta_6 >= 0, etc.
ones <- diag(1, p)
d1 <- ones[-1, ] - ones[-10,]
d1[-(1:4),] <- 0
d2 <- ones[-10,] - ones[-1,]
d2[1:4,] <- 0
d <- d1 + d2

constraints <- list(d %*% beta >= 0)
problem <- Problem(objective, constraints)
result <- solve(problem1)
result$getValue(beta)