Non positive semi-definite P matrix causes R to crash

[travis-leith]

R version 4.4.1
Osqp version 0.6.3.3
Windows 11

The following code works fine on my machine

library(conflicted)
library(tidyverse)
library(osqp)
library(Matrix)

rm(list = ls())

set.seed(123)
n <- 100
cv <- local({
  k <- 10
  
  beta <- rnorm(k * n) |> matrix(ncol = k, nrow = n)
  factor_covar <- runif(k) |> diag()
  
  idio_var <- runif(n) |> diag()
  
  beta %*% factor_covar %*% t(beta) + idio_var
})

b <- runif(n)
b <- b / sum(b)

w0 <- runif(n)
w0 <- w0 / sum(w0)

coef_b <- -runif(n)/100

aum <- 10000
b <- b * aum
w0 <- w0 * aum

P <- rbind(
  cbind(cv, matrix(0, ncol = n, nrow = n)),
  cbind(matrix(0, ncol = n, nrow = n), diag(coef_b))
) |> Matrix(sparse = TRUE)

q <- c(rep(0, 2*n))

A <- rbind(
  c(rep(1, n), rep(0, n)),
  cbind(diag(n), -diag(n)) 
) |> Matrix(sparse = TRUE)

l <- c(0, w0 - b)
u <- c(0, w0 - b)

pars <- osqpSettings(eps_abs = 1e-8, eps_rel = 1e-8)
opt <- solve_osqp(P, q, A,l, u, pars = pars)

But if I change the value of coef_b to coef_b <- -runif(n)/10, then R crashes.

Full crashing code:

library(conflicted)
library(tidyverse)
library(osqp)
library(Matrix)

rm(list = ls())

set.seed(123)
n <- 100
cv <- local({
  k <- 10
  
  beta <- rnorm(k * n) |> matrix(ncol = k, nrow = n)
  factor_covar <- runif(k) |> diag()
  
  idio_var <- runif(n) |> diag()
  
  beta %*% factor_covar %*% t(beta) + idio_var
})

b <- runif(n)
b <- b / sum(b)

w0 <- runif(n)
w0 <- w0 / sum(w0)

coef_b <- -runif(n)/10

aum <- 10000
b <- b * aum
w0 <- w0 * aum

P <- rbind(
  cbind(cv, matrix(0, ncol = n, nrow = n)),
  cbind(matrix(0, ncol = n, nrow = n), diag(coef_b))
) |> Matrix(sparse = TRUE)

q <- c(rep(0, 2*n))

A <- rbind(
  c(rep(1, n), rep(0, n)), #budget constraint
  cbind(diag(n), -diag(n)) #trade variables
) |> Matrix(sparse = TRUE)

l <- c(0, w0 - b)
u <- c(0, w0 - b)

pars <- osqpSettings(eps_abs = 1e-8, eps_rel = 1e-8)
opt <- solve_osqp(P, q, A,l, u, pars = pars)

[travis-leith]

I just tested it on a different PC and it also crashed

[bnaras]

The error message clearly points out the issue: your P matrix should be positive semi-definite.

Error in KKT matrix LDL factorization when computing the nonzero elements. The problem seems to be non-convexERROR in osqp_setup: KKT matrix factorization.
The problem seems to be non-convex.

That said, I see that it is not even strictly true for the example that worked. You should ensure it is, either by construction or by checking the eigenvalues. We might consider adding some checks for a future version but this is not a bug as I see it.

[travis-leith]

@bnaras indeed you are correct, the P matrix is not positive semi-definite. Not sure how I missed that. In any case, I do not get an error message. R just crashes. How did you get an error message? Are you on Windows 11, R 4.4.1?

This reminds me of #29, which was solved with a check on R side. But it seems like the fundamental problem of R crashing when the osqp core throws an error is not resolved.

[travis-leith]

Incidentally, setting linsys_solver = 1, when MKL is not available, also causes R to crash. I think this further illustrates my suspicion that any osqp core error causes R to crash.