Display system information.
versioninfo()
using Gadfly
# random seed
srand(2015790)
# Size of signal
n = 1024
# Sparsity (# nonzeros) in the signal
s = 20
# Number of samples (undersample by a factor of 8)
m = 128
# Generate and display the signal
x0 = zeros(n)
x0[rand(1:n, s)] = randn(s)
# Generate the random sampling matrix
A = randn(m, n) / m
# Subsample by multiplexing
y = A * x0
# plot the true signal
plot(x=1:n, y=x0, Geom.line, Guide.title("True signal x_0"))
Gurobi.jl)¶Gurobi model formulation:
objective: (1/2) x' H x + f' x
s.t. A x <= b
Aeq x = beq
lb <= x <= ub
Refer to Gurobi reference manual for setting model and parameters.
using Gurobi
env = Gurobi.Env()
setparams!(env, OutputFlag=1) # display log
# Construct the model
model = gurobi_model(env;
name = "cs",
f = ones(2 * n),
Aeq = [A -A],
beq = y,
lb = zeros(2 * n))
# Run optimization
optimize(model)
# Show results
sol = get_solution(model)
xsol = sol[1:n] - sol[n + 1:end]
plot(x=1:n, y=x0, Geom.point)
plot(x=1:n, y=xsol, Geom.line, Guide.title("Reconstructed signal overlayed with x0"))
Convex.jl¶using Convex
## Use Gurobi solver
#using Gurobi
#solver = GurobiSolver(OutputFlag=1)
#set_default_solver(solver)
# Use Mosek solver
using Mosek
solver = MosekSolver(LOG=1)
set_default_solver(solver)
## Use SCS solver
#using SCS
#solver = SCSSolver(verbose=1)
#set_default_solver(solver)
# Set up optimizaiton problem
x = Variable(n)
problem = minimize(norm(x, 1))
problem.constraints += A * x == y
# Solve the problem
@time solve!(problem)
# Display the solution
plot(x=1:n, y=x0, Geom.point)
plot(x=1:n, y=xsol, Geom.line, Guide.title("Reconstructed signal overlayed with x0"))