Skip to content

Commit 635a745

Browse files
committed
Guard vector %*% vector by whole size, not rank-2 axes
The fallthrough %*% guard (reached by vector-vector products after the gemv special cases) hardcoded rank-2 axes, emitting size(x, 2) on a rank-1 array -- a gfortran error that made conformable unknown-length dot products fail to compile. Compare whole vector sizes for rank-1 operands instead. Also document why lapack_solve()'s squareness check is routing, not a guard: rectangular solve(a, b) deliberately falls through to least squares, a tested divergence from base R.
1 parent 2c08ddd commit 635a745

3 files changed

Lines changed: 30 additions & 2 deletions

File tree

R/r2f-matrix-blas.R

Lines changed: 6 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -694,6 +694,12 @@ lapack_solve <- function(
694694

695695
nrhs <- if (b_rank == 1L) 1L else dim_or_one(B, 2L)
696696

697+
# solve(a, b) with a rectangular `a` deliberately falls through to the
698+
# least-squares branch below -- a divergence from base R (which requires
699+
# a square `a`), locked by the "least-squares" tests in
700+
# test-matrix-lapack.R. Squareness is a routing decision here, not a
701+
# correctness guard: unknown squareness routes to dgels, which solves
702+
# square systems exactly too.
697703
square <- check_conformable(m, n)
698704
if (square$ok && !square$unknown && !identical(context, "qr.solve")) {
699705
A_work <- hoist$declare_tmp(mode = "double", dims = list(m, m))

R/r2f-matrix.R

Lines changed: 6 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -108,6 +108,8 @@ register_r2f_handler(
108108
))
109109
}
110110

111+
# Vector operands (vector %*% vector reaches here) are rank-1: their
112+
# extent is their whole size, not a rank-2 axis.
111113
guard_conformable_dims(
112114
k,
113115
right_eff$rows,
@@ -116,8 +118,10 @@ register_r2f_handler(
116118
scope,
117119
left = left,
118120
right = right,
119-
left_axis = if (left_trans == "N") 2L else 1L,
120-
right_axis = if (right_trans == "N") 1L else 2L
121+
left_axis = if (left_rank == 1) NULL else if (left_trans == "N") 2L else
122+
1L,
123+
right_axis = if (right_rank == 1) NULL else if (right_trans == "N") 1L
124+
else 2L
121125
)
122126

123127
# Matrix-Matrix

tests/testthat/test-blas-guards.R

Lines changed: 18 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -43,6 +43,23 @@ test_that("triangular solve guards squareness and RHS length", {
4343
)
4444
})
4545

46+
test_that("vector %*% vector guards unknown lengths as whole sizes", {
47+
fn <- function(x, y) {
48+
declare(type(x = double(NA)), type(y = double(NA)))
49+
x %*% y
50+
}
51+
# was: the fallthrough guard hardcoded rank-2 axes, emitting size(x, 2)
52+
# on a rank-1 array -- a gfortran error that made a conformable
53+
# unknown-length dot product fail to compile at all
54+
qfn <- expect_no_warning(quick(fn))
55+
expect_equal(qfn(c(1, 2, 3), c(4, 5, 6)), c(1, 2, 3) %*% c(4, 5, 6))
56+
expect_error(
57+
qfn(c(1, 2), c(4, 5, 6)),
58+
"non-conformable arguments in %*%",
59+
fixed = TRUE
60+
)
61+
})
62+
4663
test_that("solve() guards an unknown RHS length", {
4764
fn <- function(a, b) {
4865
declare(type(a = double(2, 2)), type(b = double(NA)))
@@ -53,6 +70,7 @@ test_that("solve() guards an unknown RHS length", {
5370
expect_error(qfn(diag(2), c(1, 2, 3)), "non-conformable arguments in solve")
5471
})
5572

73+
5674
test_that("solve(a) and chol() guard squareness", {
5775
inv <- function(a) {
5876
declare(type(a = double(n, k)))

0 commit comments

Comments
 (0)