--- title: "6. Linear transformations and matrix inverse in 3D" author: "Michael Friendly" date: "`r Sys.Date()`" output: rmarkdown::html_vignette vignette: > %\VignetteIndexEntry{6. Linear transformations and matrix inverse in 3D} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} --- ```{r nomessages, echo = FALSE} knitr::opts_chunk$set( warning = FALSE, message = FALSE, echo = TRUE, collapse = TRUE, fig.height = 5, fig.width = 5 ) options(digits=4) par(mar=c(4,4,1,1)+.1) ``` ```{r setuprgl, echo=FALSE} library(rgl) library(knitr) knit_hooks$set(webgl = hook_webgl) ``` ```{r} library(matlib) # use the package library(rgl) # also using rgl ``` This vignette illustrates how linear transformations can be represented in 3D using the [`rgl` package](https://dmurdoch.github.io/rgl/). It is explained more fully in this StackExchange discussion [3D geometric interpretations of matrix inverse](https://math.stackexchange.com/questions/1948184/3d-geometric-interpretations-of-matrix-inverse). It also illustrates the determinant, $\det(\mathbf{A})$, as the volume of the transformed cube, and the relationship between $\mathbf{A}$ and $\mathbf{A}^{-1}$. ## The unit cube Start with a unit cube, which represents the identity matrix $\mathbf{I}_{3 \times 3}$ = `diag(3)` in 3 dimensions. In `rgl` this is given by `cube3d()` which returns a `"mesh3d"` object containing the vertices from -1 to 1 on each axis. I translate and scale this to make each vertex go from 0 to 1. These are represented in homogeneous coordinates to handle perspective and transformations. See `help("identityMatrix")` for details. ```{r unit-cube} 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) str(c3d) ``` The vertices are the 8 columns of the `vb` component of the object. ```{r vertices} c3d$vb ``` I want to show the transformation of $\mathbf{I}$ by a matrix $\mathbf{A}$ as the corresponding transformation of the cube. ```{r matrix-A} # matrix A: A <- matrix(c( 1, 0, 1, 0, 2, 0, 1, 0, 2), 3, 3) |> print() det(A) # A can be produced by elementary row operations on the identity matrix I <- diag( 3 ) AA <- I |> rowadd(3, 1, 1) |> # add 1 x row 3 to row 1 rowadd(1, 3, 1) |> # add 1 x row 1 to row 3 rowmult(2, 2) # multiply row 2 by 2 all(AA==A) ``` ## Define some useful functions I want to be able to display 3D mesh objects, with colored points, lines and faces. This is a little tricky, because in rgl colors are applied to vertices, not faces. The faces are colored by interpolating the colors at the vertices. See the StackOverflow discussion [drawing a cube with colored faces, vertex points and lines](https://stackoverflow.com/questions/39730889/rgl-drawing-a-cube-with-colored-faces-vertex-points-and-lines). The function `draw3d()` handles this for `"mesh3d"` objects. Another function, `vlabels()` allows for labeling vertices. ```{r draw3d} # draw a mesh3d object with vertex points and lines draw3d <- function(object, col=rep(rainbow(6), each=4), alpha=0.6, vertices=TRUE, lines=TRUE, reverse = FALSE, ...) { if(reverse) col <- rev(col) 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 $\mathbf{I}$ and $\mathbf{A}$ all together in one figure We can show the cube representing the identity matrix $\mathbf{I}$ using `draw3d()`. Transforming this by $\mathbf{A}$ is accomplished using `transform3d()`. ```{r 3d-demo1, webgl=TRUE} 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 3d-demo2, webgl=TRUE} # 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}$ The inverse of `A` is found using `solve() ```{r AI} AI <- solve(A) |> print() ``` The determinant of the matrix `A` here is 2. The determinant of $\mathbf{A}^{-1}$ is the reciprocal, $1/ \det(\mathbf{A})$. Geometrically, this means that the larger $\mathbf{A}^{-1}$ is the smaller is $\mathbf{A}^{-1}$. ```{r det-AI} det(A) det(AI) ``` You can see the relation between them by graphing both together in the same plot. It is clear that $\mathbf{A}^{-1}$ is small in the directions $\mathbf{A}$ is large, and vice versa. ```{r 3d-demo3, webgl=TRUE} 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 The `rgl` graphic can be animated by spinning it around an axis using `spin3d()`. This can be played on screen using `play3d()` or saved to a `.gif` file using `movie3d()`. ```{r 3d-demo4, webgl=TRUE} play3d(spin3d(rpm=15), duration=4) #movie3d(spin3d(rpm=15), duration=4, movie="inv-demo", dir=".") ```