Skip to content

avoid quadratic behavior when accumulating the residual when computing A*x - b #92

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.

Sign up for GitHub

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

Merged
AayushSabharwal merged 1 commit into from
Jun 4, 2026
Merged

avoid quadratic behavior when accumulating the residual when computing A*x - b #92

Changes from all commits
Commits
Show all changes
1 commit
Select commit
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
31 changes: 23 additions & 8 deletions lib/ModelingToolkitTearing/src/ModelingToolkitTearing.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -141,14 +141,7 @@ function MTKBase.unhack_system(sys::System)
additional_subs[linsolve[i]] = x[i]
end

# Compute A*x - b element-wise, skipping zero entries
for i in 1:N
res_i = -b_vec[i]
for j in 1:N
aij = A_mat[i, j]
SU._iszero(aij) && continue
res_i += aij * x[j]
end
for res_i in linear_scc_residuals(A_mat, b_vec, x)
SU._iszero(res_i) && continue
# If a linear SCC contains both `D(w)` and `w_t`, it'll contain the equation `D(w) ~ w_t`.
# When unhacking it, `D(w)` will be totermed into `w_t`. Avoid adding the `0 ~ 0` equations.
Expand Down Expand Up @@ -180,6 +173,28 @@ function MTKBase.unhack_system(sys::System)
return newsys
end

"""
$TYPEDSIGNATURES

Compute the residuals `r = A * x - b`.
"""
function linear_scc_residuals(A_mat::Matrix{SymbolicT}, b_vec::Vector{SymbolicT}, x::Vector{SymbolicT})
N = length(b_vec)
residuals = Vector{SymbolicT}(undef, N)
add_buffer = SU.ArgsT{VartypeT}()
for i in 1:N
empty!(add_buffer)
push!(add_buffer, -b_vec[i])
for j in 1:N
aij = A_mat[i, j]
SU._iszero(aij) && continue
push!(add_buffer, aij * x[j])
end
residuals[i] = SU.add_worker(VartypeT, add_buffer)
end
return residuals
end

function populate_inline_scc_map!(
inline_linear_scc_map::Dict{SymbolicT, Vector{Int}},
eq::Equation, eq_i::Int,
Expand Down
Loading