copy stuff from matlib-devel to dev/

dev/3d-demo.Rmd

@@ -0,0 +1,106 @@
+---
+title: "3D demonstrations of linear transformations and matrix inverse"
+author: "Michael Friendly"
+date: "30 Sep 2016"
+---
+Start with a unit cube, representing the identity matrix. Show its transformation
+by a matrix $A$ as the corresponding transformation of the cube.
+
+This also illustrates the determinant, det(A), as the volume of the transformed
+cube, and the relationship between $A$ and $A^{-1}$.
+
+```{r }
+library(rgl)
+library(matlib)
+# cube, with each face colored differently
+colors <- rep(2:7, each=4)
+c3d <- cube3d()
+# make it a unit cube at the origin
+c3d <- scale3d(translate3d(c3d, 1, 1, 1),
+               .5, .5, .5)
+
+# matrix A: 
+A <- matrix(c( 1, 0, 1, 0, 2, 0,  1, 0, 2), 3, 3)
+det(A)
+
+# same as the elementary row operations: 2*y, 2*z, x+y
+I <- diag( 3 )
+AA <- rowmult(rowadd(rowadd(I, 3, 1, 1), 1, 3, 1), 2, 2)
+all(AA==A)
+```
+
+## Define some useful functions
+
+```{r }
+# draw a mesh3d object with vertex points and lines
+# see: http://stackoverflow.com/questions/39730889/rgl-drawing-a-cube-with-colored-faces-vertex-points-and-lines
+draw3d <- function(object, col=rep(rainbow(6), each=4), alpha=0.6,  vertices=TRUE, lines=TRUE, ...) {
+	shade3d(object, col=col, alpha=alpha, ...)
+	vert <- t(object$vb)
+	indices <- object$ib
+	if (vertices) points3d(vert, size=5)
+	if (lines) {
+	for (i in 1:ncol(indices))
+		lines3d(vert[indices[,i],])
+		}
+}
+
+# label vertex points
+vlabels <- function(object, vertices, labels=vertices, ...) {
+	text3d( t(object$vb[1:3, vertices] * 1.05), texts=labels, ...)
+}
+```
+
+## show I and A all together in one figure
+
+```{r }
+open3d()
+draw3d(c3d)
+vlabels(c3d, c(1,2,3,5))
+
+axes <- rbind( diag(3), -diag(3) )
+rownames(axes) <- c("x", "y", "z", rep(" ", 3))
+vectors3d(axes, frac.lab=1.2, headlength = 0.2, radius=1/20, lwd=3)
+
+c3t<- transform3d(c3d, A) 
+draw3d(c3t)
+vlabels(c3t, c(1,2,3,5))
+```
+
+## Same, but using separate figures, shown side by side
+
+```{r }
+# NB: this scales each one separately, so can't see relative size
+open3d()
+mfrow3d(1,2, sharedMouse=TRUE)
+draw3d(c3d)
+vectors3d(axes, frac.lab=1.2, headlength = 0.2, radius=1/20, lwd=3)
+
+next3d()	
+draw3d(c3t)
+vectors3d(axes, frac.lab=1.2, headlength = 0.2, radius=1/20, lwd=3)
+```
+
+## A and A^{-1}
+
+```{r }
+open3d()
+draw3d(c3t)
+c3Inv <- transform3d(c3d, solve(A))
+draw3d(c3Inv)
+vectors3d(axes, frac.lab=1.2, headlength = 0.2, radius=1/20, lwd=3)
+vlabels(c3t, 8, "A", cex=1.5)
+vlabels(c3t, c(2,3,5))
+
+vlabels(c3Inv, 4, "Inv", cex=1.5)
+vlabels(c3Inv, c(2,3,5))
+```
+
+Animate
+
+```{r }
+play3d(spin3d(rpm=15),  duration=4)
+
+movie3d(spin3d(rpm=15), duration=4, movie="inv-demo", dir=".")
+```
+

dev/SVD-demo.Rmd

@@ -0,0 +1,82 @@
+---
+title: "Visualizing the SVD"
+author: "Duncan Murdoch"
+date: "February 2, 2015"
+output: 
+  html_document
+---
+
+These notes are taken from a presentation by Duncan Murdoch, titled "Recent Developments in rgl".
+
+## The Singular Value Decomposition
+
+For a square $n \times n$ matrix $A$, the SVD is
+$$ A = U D V^T $$
+where 
+
+- $U$ and $V$ are $n \times n$ orthogonal matrices (i.e. $U^TU = V^TV = I$)
+- $D$ is an $n \times n$ diagonal matrix with non-negative
+entries
+
+## Displaying a Matrix Graphically
+
+- Matrices are representations of
+linear operators on vector spaces.
+- The matrix $A$ is characterized by the behaviour of $y = A x$
+as we vary $x$.
+- Use the `rgl` package to develop a graphical representation of 
+$3 \times 3$ matrices, e.g.
+$$ A = \left(\begin{array}{rrr}
+1 & 0.1 & 0.1 \\
+2 & 1 & 0.1 \\
+0.1 & 0.1 & 0.5
+\end{array}\right) $$
+
+## Attempt 1:  display the actions on the basis vectors
+
+```{r rgl=TRUE,fig.height=2}
+A <- matrix(c(1,2,0.1, 0.1,1,0.1, 0.1,0.1,0.5), 3,3)
+basis <- cbind(c(0,1,0,0,0,0), c(0,0,0,1,0,0), c(0,0,0,0,0,1))
+segments3d(basis, lwd = 3)
+segments3d(basis %*% t(A), col = "red", lwd = 5)
+text3d(1.1*basis, texts = c("","x", "", "y", "", "z"), cex = 2)
+text3d(1.1*basis %*% t(A), col = "red", 
+       texts=c("", "x", "", "y", "", "z"), cex = 2)
+```
+
+## Attempt 2:  coloured ellipsoids.
+
+```{r rgl=TRUE, echo=FALSE}
+sphere <- subdivision3d(cube3d(color=rep(rainbow(6),rep(4*4^4,6))),depth=4)
+sphere$vb[4,] <- apply(sphere$vb[1:3,], 2, function(x) sqrt(sum(x^2)))
+mult <- function(matrix, obj) transform3d(obj, t(matrix))
+mat <- matrix(1:6, 2, 3)
+layout3d(rbind(mat,mat+6), height = c(3,1,3,1), sharedMouse=TRUE)
+shade3d(sphere)
+next3d()
+text3d(0,0,0, "Identity", cex=1.5)
+next3d()
+next3d(reuse=FALSE)
+next3d(reuse=FALSE)
+shade3d(mult(A, sphere))
+next3d()
+text3d(0,0,0,"A", cex=1.5)
+
+svd <- svd(A)
+U <- svd$u
+D <- diag(svd$d)
+V <- svd$v
+next3d()
+shade3d(mult(U, sphere))
+next3d()
+text3d(0,0,0, "U", cex=1.5)
+next3d()
+shade3d(mult(D, sphere))
+next3d()
+text3d(0,0,0, "D", cex=1.5)
+next3d()
+shade3d(mult(V, sphere))
+next3d()
+text3d(0,0,0, "V", cex=1.5)
+```
+