base::chol() slows down when matrix contains many small entries

Question

I've noticed that base::chol() severely slows down when the matrix contains many small elements. Here is an example:

## disable openMP
library(RhpcBLASctl); blas_set_num_threads(1); omp_set_num_threads(1)

• Baseline: create positive definite matrix and get timing for chol().

  loc <- expand.grid(1:60, 1:50)
  covmat1 <- exp(-as.matrix(dist(loc)))
  mean(c(covmat1))
  # [1] 0.002076862
  system.time(chol1 <- chol(covmat1))
  #   user  system elapsed 
  #  0.313   0.024   0.337 
  

• Increase small values: create covmat2 matrix with more small values.

  covmat2 <- exp(-as.matrix(dist(loc))*10)
  mean(c(covmat2))
  # [1] 0.0003333937
  system.time(chol2 <- chol(covmat2))
  #   user  system elapsed 
  #  2.311   0.021   2.333 
  

  Compared to the base line this slows down the computation by almost factor 10.

• Set small values to zero: set values of covmat2 that are smaller than 1e-13 to zero.

  covmat3 <- covmat2
  covmat3[covmat3 < 1e-13] <- 0
  mean(c(covmat3))
  # [1] 0.0003333937
  system.time(chol3 <- chol(covmat3))
  #   user  system elapsed 
  #  0.302   0.016   0.318 
  

  This version is again faster and similar to the base line.

Why does this slowdown happen?

---

Notes:

Repeated evaluations of the timing experiments lead to similar results.

I know that for matrices with many values close to zero it might be more efficient to use a sparse matrix approach, e.g., the R package spam.

sessionInfo()
## R version 3.6.1 (2019-07-05)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Linux Mint 19.2

## Matrix products: default
## BLAS:   /usr/lib/x86_64-linux-gnu/openblas/libblas.so.3
## LAPACK: /usr/lib/x86_64-linux-gnu/libopenblasp-r0.2.20.so

Armoury

Strange but maybe the running time increased because of the numbers. You can check also Rfast's cholesky but they use omp for parallel computing. – Immense

@Csd in my understanding the computation time of the Cholesky should only depend on the dimension of the matrix and not on the numbers in it. – Armoury

What is the running time for multiplying two numbers? n*m where n,m are the number of digits of each number. It is although strange but it could be explained in a way. Some systems when there is a multiplication that is zero they stop. It is an optimization that I have learn in school. That is why when there are many zeros it is fastest than before. – Immense

?chol the third line under details suggests that there may be some early checks to determine whether the values are zero or negative. Just a guess without looking at the source code – Stilt

@Csd Are you sure that the execution time to multiply two variables of type double depends on what number is stored? This would be surprising to me. But I am happy to look at a reference if you have one. - Even if this would be true, it is not likely that the binary representations of the numbers in covmat1 is shorter than those of covmat2. Hence, I don't think this could explain the timing difference. – Armoury

As I told you, it might be your cpu that optimizes the calculations. In university we have learned that some CPUs do that. For example if you want to multiply two numbers and one of them is zero your software will return immediately zero. A multiplication of two numbers that one of them is zero leads to zero result but a multiplication of two numbers to range 1e-13 will lead to some number. This means the cpu multiply indeed the numbers. If a multiplication of type double has runtime x then count the total x you will need for all the multiplications. – Immense

@Csd I see what you mean. Thanks for the comment. But covmat2 has more small numbers than covmat1 and according to your explanation it should be faster. But it is slower... – Armoury

According to my answer covmat3 is faster which is indeed... – Immense

Small number in memory doesn't mean zero. If you multiplying two small numbers you will get a small number but if you multiply zero with a small number you will get zero. – Immense

@ManosPapadakis I see. I think this is the correct explanation. – Armoury