Are the matrix operations in R programming optimized for the matrices of specific shapes (e.g. banded matrices)?

• Theoretically, certain matrix operations (e.g., multiplication, inverse, eigendecomposition) are more efficient when they are applied to some special matrices (e.g. banded matrix). Hence, I wonder whether the base-R matrix operators (%*%) and functions (solve, eigen) utilizes the structure of specific matrices for faster computing. It seems the answer is no, see the following example:

set.seed(2501)
A <- B <- matrix(runif(1000*1000), 1000, 1000)
for (i in 1:1000) {
  idx <- (i-2):(i+2)
  idx <- idx[idx >= 1 & idx <= 1000]
  A[i,-idx] <- 0
}

X <- matrix(runif(1000*1000), 1000, 1000)
system.time(X %*% A %*% t(X))
#> user system elapse 
#> 2.25 0.03 2.30
system.time(X %*% B %*% t(X))
#> user system elapse 
#> 2.24 0.03 2.28

^{Created on 2022-06-21 by the reprex package (v2.0.1)}

• If not, are there any user-friendly APIs that are optimized for matrix operations?