I am trying to solve a MIQP portfolio optimisation problem using JuMP and SCIP (tried Juniper before but it crashed Julia). There are repeated errors:
Error: no BLAS/LAPACK library loaded for dsyev_()
[nlpi_ipopt.cpp:2680] ERROR: There was an error when calling DSYEV. INFO = 3
# Afterwards
WARNING: Failed to compute eigenvalues of quadratic coefficient matrix --> don't know curvature
which I believe may imply that SCIP’s IPOPT plugin can’t find or interface with LAPACK at runtime — even though Julia itself can. Indeed we can see lbt as the result of
LinearAlgebra.BLAS.vendor()
Julia is using the Apple Accelerate framework on macOS. So in principle it should support LAPACK routines like dsyev_(). I can check curvature myself with something like
Sigma = Symmetric(Sigma)
Sigma += 1e-6 * I # regularization to ensure positive definiteness
eigvals = eigen(Sigma).values
minimum(eigvals)
Is there a way to tell JuMP/SCIP what the eigenvalues are since, as it is, it is unable to compute them? Or perhaps there is a better way to overcome this problem? Thank you!
Thank you very much for having a look at this. This is the exact script, not very minimal, but reproduces the error. Notice that the script runs till the end. I uploaded the CSV file with the DJIA returns (added the .jl extension (djia.csv.jl) so that I could upload the file …).
############## Cardinality Constraints ######################
using LinearAlgebra
using JuMP
using SCIP
using CSV
# Get returns from Dow Jones 30 components from 2024-01-01, to 2025-03-31
returns = CSV.read("djia.csv", DataFrame) # Load the data
# Calculate sample mean vector (mu) and covariance matrix (Sigma)
mu = mean(Matrix(returns[:, 2:end]), dims=1)'
n = length(mu)
Sigma = cov(Matrix(returns[:, 2:end]))
# Market capitalization weights (assumed equal for simplicity)
w_mkt = fill(1/n, n)
# Risk aversion coefficient
lambda = 5
# Implied equilibrium returns
PI = lambda * Sigma * w_mkt
# Investor views
# Example: Expect AAPL to outperform MSFT by 2%
P = [1 -1 zeros(1, n-2)]
q = [0.02]
# Uncertainty of views
Omega = P * Sigma * P' * 0.25
# Black-Litterman formula
tau = 0.05
M = inv(inv(tau * Sigma) + P' * inv(Omega) * P)
mu_bl = M * (inv(tau * Sigma) * PI + P' * inv(Omega) * q)
# Create the SCIP/JuMP model
model = Model(SCIP.Optimizer)
# Weights: long-only, max 10%
@variable(model, 0 <= w[1:n] <= 0.1)
# Binary variables for activation
@variable(model, z[1:n], Bin)
# Fully invested portfolio
@constraint(model, sum(w) == 1.0)
# Threshold if invested
@constraint(model, [i=1:n], w[i] >= 0.01 * z[i])
# Upper bound when active
@constraint(model, [i=1:n], w[i] <= 0.1 * z[i])
# Objective: use Black-Litterman (mu_bl) or other expected return vector
objective = dot(mu_bl, w) - (lambda/2) * dot(w, Sigma * w)
@objective(model, Max, objective)
optimize!(model)
w_dc = value.(w)
println("Optimal weights (Disjunctive Constraint) = ", round.(w_dc,digits = 6))
I can’t reproduce this, so try updating your packages. There have been a bunch of improvements between SCIP 0.11.6 (released in October 2022) and the latest version. It’s probably something to do with your M1 Mac.
Indeed, updating solved the problem. I was not aware I had such an old version installed. Also tried with Juniper + OSQP and it worked fine too. Thank you!