What happens
solve(a, b) compiles for rectangular a and returns the least-squares
solution — a value where R raises an error:
fn <- function(X, y) {
declare(type(X = double(n, k)), type(y = double(n)))
solve(X, y)
}
set.seed(1)
X <- matrix(rnorm(6), 3, 2)
y <- rnorm(3)
solve(X, y)
#> Error in solve.default(X, y) : 'a' (3 x 2) must be square
quick(fn)(X, y)
#> [1] -0.5511845 0.1229073 # qr.solve(X, y), silently
Why it happens
lapack_solve() in R/r2f-matrix-blas.R treats squareness as a routing
decision, not a requirement: statically square systems go to LU (dgesv),
everything else — statically rectangular and unknown-at-compile-time —
falls through to a least-squares dgels call. So a symbolic-dims solve()
silently behaves like qr.solve() whenever the runtime matrix happens to
be rectangular.
Expected behavior
R's solve() requires a square a; least squares is qr.solve()'s job,
and quickr already supports qr.solve() separately. A statically
rectangular a should be a compile error, and unknown squareness should
be checked at run time before the dgesv call — the same
compile-error/runtime-guard policy the other linear-algebra lowerings
follow.
What happens
solve(a, b)compiles for rectangularaand returns the least-squaressolution — a value where R raises an error:
Why it happens
lapack_solve()in R/r2f-matrix-blas.R treats squareness as a routingdecision, not a requirement: statically square systems go to LU (dgesv),
everything else — statically rectangular and unknown-at-compile-time —
falls through to a least-squares dgels call. So a symbolic-dims
solve()silently behaves like
qr.solve()whenever the runtime matrix happens tobe rectangular.
Expected behavior
R's
solve()requires a squarea; least squares isqr.solve()'s job,and quickr already supports
qr.solve()separately. A staticallyrectangular
ashould be a compile error, and unknown squareness shouldbe checked at run time before the dgesv call — the same
compile-error/runtime-guard policy the other linear-algebra lowerings
follow.