R Under development (unstable) (2024-12-10 r87437 ucrt) -- "Unsuffered Consequences"
Copyright (C) 2024 The R Foundation for Statistical Computing
Platform: x86_64-w64-mingw32/x64

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> #### NOT part of the cluster package!
> 
> #### Find out what exactly sweep() in ../src/spannel.f is doing
> #### in order to eventually replace it with BLAS calls !
> 
> ### subroutine sweep (cov,nord,ixlo,nel,deter)
> ###            ===============================
> ### is called only once as
> ###	call sweep(cov,ndep,0,i,deter)
> ### where  i in 0:ndep
> 
> sweep1 <- function(cov, i, det = 1)
+ {
+   ## Purpose:
+   ## -------------------------------------------------------------------------
+   ## Arguments:
+   ## -------------------------------------------------------------------------
+   ## Author: Martin Maechler, Date: 22 Jan 2002, 08:58
+     if(!is.matrix(cov) || 0 != diff(D <- dim(cov)))
+         stop("'cov' must be a square matrix")
+     if((nord <- as.integer(D[1] - 1)) < 1)## cov[0:nord, 0:nord]
+         stop("'cov' must be at least 2 x 2")
+     if(0 > (i <- as.integer(i)) || i > nord)
+         stop("'i' must be in 0:nord, where nord = nrow(cov)-1")
+     storage.mode(cov) <- "double"
+     .C(cluster:::cl_sweep,
+        cov,
+        nord,
+        ixlo = 0:0,
+        i = i,
+        deter=det)
+ }
> 
> sweepAll <- function(cov, det = 1)
+ {
+   ## Purpose:
+   ## -------------------------------------------------------------------------
+   ## Arguments:
+   ## -------------------------------------------------------------------------
+   ## Author: Martin Maechler, Date: 22 Jan 2002, 08:58
+     if(!is.matrix(cov) || 0 != diff(D <- dim(cov)))
+         stop("'cov' must be a square matrix")
+     if((nord <- as.integer(D[1] - 1)) < 1)## cov[0:nord, 0:nord]
+         stop("'cov' must be at least 2 x 2")
+     storage.mode(cov) <- "double"
+     for(i in 0:nord) {
+ 	.C(cluster:::cl_sweep,
+ 	   cov,
+ 	   nord,
+ 	   ixlo = 0:0,
+ 	   i = i,
+ 	   deter = det,
+ 	   DUP = FALSE) # i.e. work on 'cov' and 'det' directly
+         if(det <= 0)
+             cat("** i = ", i, "; deter = ", format(det)," <= 0\n",sep="")
+     }
+     list(cov = cov, deter = det)
+ }
> 
> require(cluster)
Loading required package: cluster
> 
> ## Examples with errors
> m1 <- cov(cbind(1, 1:5))
> try(sweepAll(m1))# deter = 0; cov[2,2] = Inf
$cov
     [,1] [,2]
[1,]    0  0.0
[2,]    0  2.5

$deter
[1] 1

> 
> ## ok
> (m2 <- cov(cbind(1:5, c(2:5,1), c(4:2,2,6))))
     [,1] [,2]  [,3]
[1,] 2.50  0.0  0.75
[2,] 0.00  2.5 -2.50
[3,] 0.75 -2.5  2.80
> qr(m2, tol = .001)$rank
[1] 3
> sweepAll(m2) ## deter = 0
$cov
     [,1] [,2]  [,3]
[1,] 2.50  0.0  0.75
[2,] 0.00  2.5 -2.50
[3,] 0.75 -2.5  2.80

$deter
[1] 1

> 
> proc.time()
   user  system elapsed 
   0.21    0.12    0.28