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

Update semiquadratic_form to use more intuitive form #1131

Open
4 tasks
bowenszhu opened this issue May 6, 2024 · 0 comments
Open
4 tasks

Update semiquadratic_form to use more intuitive form #1131

bowenszhu opened this issue May 6, 2024 · 0 comments

Comments

@bowenszhu
Copy link
Member

bowenszhu commented May 6, 2024

"""
$(TYPEDSIGNATURES)
Returns a tuple of 4 objects:
1. a matrix `A` of dimensions (m x n)
2. a matrix `B` of dimensions (m x (n+1)*n/2)
3. a vector `v2` of length (n+1)*n/2 containing monomials of `vars` upto degree 2 and zero where they are not required.
4. a residual vector `c` of length m.
where `n == length(exprs)` and `m == length(vars)`.
The result is arranged such that, `A * vars + B * v2 + c` is the same as `exprs`.
"""
function semiquadratic_form(exprs, vars)
ds, nls = semipolynomial_form(exprs, vars, 2; consts = false)
idxmap = Dict(v=>i for (i, v) in enumerate(vars))
m, n = length(exprs), length(vars)
I1 = Int[]
J1 = Int[]
V1 = Num[]
I2 = Int[]
J2 = Int[]
V2 = Num[]
v2_I = Int[]
v2_V = Num[]
for (i, d) in enumerate(ds)
for (k, v) in d
if pdegree(k) == 1
push!(I1, i)
push!(J1, idxmap[k])
push!(V1, v)
elseif pdegree(k) == 2
push!(I2, i)
if isop(k, ^)
b, e = arguments(k)
@assert e == 2
q = idxmap[b]
j = div(q*(q+1), 2)
push!(J2, j) # or div(q*(q-1), 2) + q
push!(V2, v)
else
@assert isop(k, *)
a, b = unsorted_arguments(k)
p, q = extrema((idxmap[a], idxmap[b]))
j = div(q*(q-1), 2) + p
push!(J2, j)
push!(V2, v)
end
push!(v2_I, j)
push!(v2_V, k)
else
error("This should never happen")
end
end
end
#v2 = SparseVector(div(n * (n + 1), 2), v2_I, v2_V) # When it works in the future
# until then
v2 = zeros(Num, div(n * (n + 1), 2))
v2[v2_I] .= v2_V
tuple(sparse(I1,J1,V1, m, n),
sparse(I2,J2,V2, m, div(n * (n + 1), 2)),
v2,
wrap.(nls))
end

Description

Currently, semiquadratic_form returns a tuple of four objects: A, B, v2, and c arranged such that A * vars + B * v2 + c is equivalent to the input expression.
This issue proposes a breaking change to update the function's behavior to a more intuitive and user-friendly form:
A * vars + vars' * B * vars + c

Motivation

The current form, while functional, requires constructing a vector v2 containing monomials of vars up to degree 2. This approach is less intuitive and can be cumbersome to work with. The proposed change aligns with standard mathematical notation for quadratic forms and simplifies subsequent calculations.

Proposed Changes

  • Modify the internal logic to directly compute matrices B in the new form.
  • Update the function's docstring.
  • Adjust existing tests to accommodate the new output format.
  • Provide additional examples demonstrating the usage of the updated function in documentation.

Impact

This change is considered breaking as it modifies the output structure of semiquadratic_form. However, considering the function's limited usage, the impact on existing code is expected to be minimal.

cc @shashi

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

1 participant