Skip to content

Instantly share code, notes, and snippets.

@friendly
Created Sep 30, 2016
Embed
What would you like to do?
3D demonstrations of linear transformations and matrix inverse
#' ---
#' 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}$.
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
# 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
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
# 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}
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
play3d(spin3d(rpm=15), duration=4)
movie3d(spin3d(rpm=15), duration=4, movie="inv-demo", dir=".")
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment