Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Move isfeasible to util section #3759

Open
wants to merge 1 commit into
base: master
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 0 additions & 1 deletion src/Interfaces/AbstractPolyhedron.jl
Original file line number Diff line number Diff line change
Expand Up @@ -76,7 +76,6 @@ function _linear_map_hrep_helper(M::AbstractMatrix, P::LazySet,
end

# internal functions; defined here due to dependency SymEngine and submodules
function isfeasible end
function remove_redundant_constraints end
function remove_redundant_constraints! end
function _ishalfspace end
Expand Down
44 changes: 0 additions & 44 deletions src/Interfaces/AbstractPolyhedron_functions.jl
Original file line number Diff line number Diff line change
@@ -1,7 +1,6 @@
export constrained_dimensions,
remove_redundant_constraints,
remove_redundant_constraints!,
isfeasible,
addconstraint!,
ishyperplanar

Expand Down Expand Up @@ -952,49 +951,6 @@ function _project_polyhedron(P::LazySet, block; kwargs...)
return constraints_list(πP)
end

"""
isfeasible(A::AbstractMatrix, b::AbstractVector, [witness]::Bool=false;
[solver]=nothing)

Check for feasibility of linear constraints given in matrix-vector form.

### Input

- `A` -- constraints matrix
- `b` -- constraints vector
- `witness` -- (optional; default: `false`) flag for witness production
- `solver` -- (optional; default: `nothing`) LP solver

### Output

If `witness` is `false`, the result is a `Bool`.

If `witness` is `true`, the result is a pair `(res, w)` where `res` is a `Bool`
and `w` is a witness point/vector.

### Algorithm

This implementation solves the corresponding feasibility linear program.
"""
function isfeasible(A::AbstractMatrix, b::AbstractVector, witness::Bool=false;
solver=nothing)
N = promote_type(eltype(A), eltype(b))
# feasibility LP
lbounds, ubounds = -Inf, Inf
sense = '<'
obj = zeros(N, size(A, 2))
if isnothing(solver)
solver = default_lp_solver(N)
end
lp = linprog(obj, A, sense, b, lbounds, ubounds, solver)
if is_lp_optimal(lp.status)
return witness ? (true, lp.sol) : true
elseif is_lp_infeasible(lp.status)
return witness ? (false, N[]) : false
end
return error("LP returned status $(lp.status) unexpectedly")
end

# convenience function to invert the result of `isfeasible` while still including the witness result
function _isinfeasible(constraints::AbstractVector{<:HalfSpace},
witness::Bool=false; solver=nothing)
Expand Down
2 changes: 0 additions & 2 deletions src/Sets/HalfSpace/isfeasible.jl
Original file line number Diff line number Diff line change
@@ -1,5 +1,3 @@


"""
isfeasible(constraints::AbstractVector{<:HalfSpace}, [witness]::Bool=false;
[solver]=nothing)
Expand Down
45 changes: 45 additions & 0 deletions src/Utils/lp_solvers.jl
Original file line number Diff line number Diff line change
@@ -1,3 +1,5 @@
export isfeasible

function default_lp_solver(::Type{T}) where {T}
key = task_local_lp_solver_key(T)
LP = get!(() -> JuMP.Model(default_lp_solver_factory(T)), task_local_storage(), key)
Expand Down Expand Up @@ -35,3 +37,46 @@ function default_lp_solver_polyhedra(N; kwargs...)
require(@__MODULE__, :Polyhedra; fun_name="default_lp_solver_polyhedra")
return error("no default solver for numeric type $N")
end

"""
isfeasible(A::AbstractMatrix, b::AbstractVector, [witness]::Bool=false;
[solver]=nothing)

Check for feasibility of linear constraints given in matrix-vector form.

### Input

- `A` -- constraints matrix
- `b` -- constraints vector
- `witness` -- (optional; default: `false`) flag for witness production
- `solver` -- (optional; default: `nothing`) LP solver

### Output

If `witness` is `false`, the result is a `Bool`.

If `witness` is `true`, the result is a pair `(res, w)` where `res` is a `Bool`
and `w` is a witness point/vector.

### Algorithm

This implementation solves the corresponding feasibility linear program.
"""
function isfeasible(A::AbstractMatrix, b::AbstractVector, witness::Bool=false;
solver=nothing)
N = promote_type(eltype(A), eltype(b))
# feasibility LP
lbounds, ubounds = -Inf, Inf
sense = '<'
obj = zeros(N, size(A, 2))
if isnothing(solver)
solver = default_lp_solver(N)
end
lp = linprog(obj, A, sense, b, lbounds, ubounds, solver)
if is_lp_optimal(lp.status)
return witness ? (true, lp.sol) : true
elseif is_lp_infeasible(lp.status)
return witness ? (false, N[]) : false
end
return error("LP returned status $(lp.status) unexpectedly")
end
Loading