Skip to main content

powerio_matrix/matrix/
lacpf.rs

1//! Linear AC power flow (LACPF) block matrix.
2//!
3//! At flat start `u* = 1 + j0`, the linearized power flow Jacobian is
4//!
5//! ```text
6//!   F(x*) = [[ G  -B  -I  0 ],
7//!            [-B  -G   0 -I ]],
8//! ```
9//!
10//! whose `2n × 2n` block (without the load injection identity columns) is
11//!
12//! ```text
13//!   J = [[ G  -B ],
14//!        [-B  -G ]].
15//! ```
16//!
17//! Satisfies `p = G ε - B θ`, `q = -B ε - G θ`. Indefinite (saddle point);
18//! emitted as a hard input alongside Bp, Bpp, and ±Im(Y_bus).
19
20use sprs::CsMat;
21
22use crate::Result;
23use crate::indexed::IndexedNetwork;
24
25use super::BuildOptions;
26use super::ybus::calc_admittance_matrix;
27
28/// Calculate the flat start linear AC power flow block matrix.
29pub fn calc_lacpf_matrix(case: &IndexedNetwork, opts: &BuildOptions) -> Result<CsMat<f64>> {
30    let parts = calc_admittance_matrix(case, opts)?;
31    let n = case.n();
32    let two_n = 2 * n;
33
34    // Walk both G and B once and emit the 4 blocks: [+G, -B; -B, -G].
35    // `calc_admittance_matrix` already returns CSR, so iterate the parts
36    // directly rather than deep-copying through `to_csr()`.
37    let mut tri = sprs::TriMat::with_capacity((two_n, two_n), 2 * (parts.g.nnz() + parts.b.nnz()));
38
39    for (&v, (i, j)) in &parts.g {
40        tri.add_triplet(i, j, v); // top-left:  +G
41        tri.add_triplet(n + i, n + j, -v); // bottom-right: -G
42    }
43    for (&v, (i, j)) in &parts.b {
44        tri.add_triplet(i, n + j, -v); // top-right:    -B
45        tri.add_triplet(n + i, j, -v); // bottom-left:  -B
46    }
47
48    Ok(tri.to_csr())
49}