Skip to content

Commit

Permalink
Merge pull request #6 from darnstrom/python
Browse files Browse the repository at this point in the history
Plotting + x->z
  • Loading branch information
darnstrom authored Nov 5, 2024
2 parents 9e5b81f + d2300fd commit e5239da
Show file tree
Hide file tree
Showing 6 changed files with 39 additions and 19 deletions.
6 changes: 3 additions & 3 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -4,15 +4,15 @@

$$
\begin{align}
\min_{x} & ~\frac{1}{2}x^{T}Hx+(f+F \theta)^{T}x \\
\text{s.t.} & ~A x \leq b + B \theta \\
\min_{z} & ~\frac{1}{2}z^{T}Hz+(f+F \theta)^{T}x \\
\text{s.t.} & ~A z \leq b + B \theta \\
& ~\theta \in \Theta
\end{align}
$$

where $H \succ 0$ and $\Theta \triangleq \lbrace l \leq \theta \leq u : A_{\theta} \theta \leq b_{\theta}\rbrace$.

The solution $x^*(\theta)$ is a piecewise-affine function over a polyhedral partition.
The solution $z^*(\theta)$ is a piecewise-affine function over a polyhedral partition.

## Example
The following code solves the mpQP in Section 7.1 in Bemporad et al. 2002
Expand Down
34 changes: 27 additions & 7 deletions src/plot.jl
Original file line number Diff line number Diff line change
Expand Up @@ -7,7 +7,7 @@ using RecipesBase
vs=PolyDAQP.vrep_2d(minrep(slice(r.Ath,r.bth,collect(3:nth))...)...)
nv = length(vs)
x,y = first.(vs),last.(vs)
z = [r.x[:,1]'*[v;zeros(nth-2);1] for v in vs]
z = [r.z[:,1]'*[v;zeros(nth-2);1] for v in vs]
connections--> (zeros(Int, nv-2),collect(1:nv-2),collect(2:nv-1))
return x,y,z
end
Expand All @@ -27,8 +27,8 @@ end
nv = length(vs)
nv < 2 && continue
x,y = first.(vs), last.(vs)
c = r.x[ids,uid]'*values + r.x[end,uid]
z = [r.x[free_ids,uid]'*v + c for v in vs]
c = r.z[ids,uid]'*values + r.z[end,uid]
z = [r.z[free_ids,uid]'*v + c for v in vs]
@series begin
st --> :mesh3d
legend --> false
Expand All @@ -41,7 +41,13 @@ end
end

## PGFPlotsX
function pplot(rs::Vector{CriticalRegion};uid=0, fix_ids = nothing, fix_vals=nothing,opts=Dict{Symbol,Any}())
function plot_regions(sol::Solution;fix_ids=nothing,fix_vals=nothing,opts=Dict{Symbol,Any}())
pplot(get_critical_regions(sol);z_id=0,fix_ids,fix_vals,opts)
end
function plot_solution(sol::Solution;z_id=1,fix_ids=nothing,fix_vals=nothing,opts=Dict{Symbol,Any}())
pplot(get_critical_regions(sol);z_id,fix_ids,fix_vals,opts)
end
function pplot(rs::Vector{CriticalRegion};z_id=0, fix_ids = nothing, fix_vals=nothing,opts=Dict{Symbol,Any}())
isempty(rs) && error("Cannot plot empty collection")
nth = size(rs[1].Ath,1)
ids = isnothing(fix_ids) ? collect(3:nth) : fix_ids
Expand All @@ -55,9 +61,23 @@ function pplot(rs::Vector{CriticalRegion};uid=0, fix_ids = nothing, fix_vals=not
p = Polyhedron(slice(r.Ath,r.bth,ids;values)...)
isempty(p) && continue
push!(ps,p)
uid == 0 && continue
c = r.x[ids,uid]'*values + r.x[end,uid]
push!(fs,v->c+r.x[free_ids,uid]'*v[1:2])
z_id == 0 && continue
c = r.z[ids,z_id]'*values + r.z[end,z_id]
push!(fs,v->c+r.z[free_ids,z_id]'*v[1:2])
end
# Some plotting
lopts = Dict(
:xlabel=>"\\large\$\\theta_{"*string(free_ids[1])*"}\$",
:ylabel=>"\\large\$\\theta_{"*string(free_ids[2])*"}\$",
)
if z_id != 0
push!(lopts,:zlabel=>"\\large\$z^*_{"*string(z_id[1])*"}\$")

lopts = merge(Dict(:view=>(45,45),),lopts)
else
push!(lopts,:ylabel_style=> "{yshift={-5pt}}")
end

opts = merge(lopts,opts)
PolyDAQP.pplot(ps;fs,opts)
end
12 changes: 6 additions & 6 deletions src/tree.jl
Original file line number Diff line number Diff line change
Expand Up @@ -47,22 +47,22 @@ end
function get_feedbacks(CRs; tol=1e-6)
isempty(CRs) && return nothing
nth = size(CRs[1].Ath,1)
X,ids = [], []
Z,ids = [], []
for cr in CRs
id = 0
for (i,x) in enumerate(X)
if norm(cr.x-x) < tol
for (i,z) in enumerate(Z)
if norm(cr.z-z) < tol
id = i
push!(ids,id)
break
end
end
if id == 0 # No duplicate
push!(X,cr.x)
push!(ids,length(X))
push!(Z,cr.z)
push!(ids,length(Z))
end
end
return X,ids
return Z,ids
end

# TODO: Can be cut in half by using points in CR
Expand Down
2 changes: 1 addition & 1 deletion src/types.jl
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@ mutable struct CriticalRegion
AS::Vector{Int16}
Ath::Matrix{Float64}
bth::Vector{Float64}
x::Matrix{Float64}
z::Matrix{Float64}
lam::Matrix{Float64}
th::Vector{Float64}
end
Expand Down
2 changes: 1 addition & 1 deletion src/utils.jl
Original file line number Diff line number Diff line change
Expand Up @@ -88,7 +88,7 @@ function denormalize(cr::CriticalRegion,scaling,translation)
else
An,bn = zeros(0,0),zeros(0)
end
xn = !isempty(cr.x) ? denormalize(cr.x,scaling,translation) : zeros(0,0)
xn = !isempty(cr.z) ? denormalize(cr.z,scaling,translation) : zeros(0,0)
lamn = !isempty(cr.lam) ? denormalize(cr.lam,scaling,translation) : zeros(0,0)
thn = !isempty(cr.th) ? denormalize(cr.th,scaling,translation) : zeros(0)
return CriticalRegion(cr.AS,An,bn,xn,lamn,thn)
Expand Down
2 changes: 1 addition & 1 deletion test/runtests.jl
Original file line number Diff line number Diff line change
Expand Up @@ -57,7 +57,7 @@ end
inds = pointlocation(θ,sol.CRs);
containment_inds[n]=length(inds)
AS = sol.CRs[inds[1]].AS
xsol = sol.CRs[inds[1]].x'*[θ;1]
xsol = sol.CRs[inds[1]].z'*[θ;1]
λsol = sol.CRs[inds[1]].lam'*[θ;1]
f = mpQP.f[:,1]+mpQP.F*θ
b = mpQP.b[:,1]+mpQP.B*θ
Expand Down

0 comments on commit e5239da

Please sign in to comment.