【R言語で作成するアニメーション】楕円・斜め楕円のアニメーションを作図してみた件

はじめに
楕円の作図
斜め楕円の作図
後書き
参考資料

はじめに

ふと、夏季休暇中に、楕円の公式のパラメータを変えたときに、どんな図形になるんだったっけと思い立ち、 Rで、アニメーションのプログラムを作ってみた。

楕円の公式を解析的に解いてみて作図したところ、楕円の形って結構発散するんだなとはじめて知った、夏休み。

楕円の作図

#インストール
install.packages("plotrix")

#楕円の公式(R風)
#x^2/a^2 + y^2/b^2 = 1
#y = ±sqrt((1 - x^2/a^2 )*b^2)

#パラメータ設定
A <- c(1:10, 10:1, 1:3)
B <- rep(3, length(A))
a <- c(A, 3, B)
b <- c(B, 2, A)

animation::saveGIF({
for(n in seq_len(length(a))){
#n <- 1
#X軸
x <- seq(-10, 10, by=0.01)

#Y軸
suppressWarnings(y1 <- sqrt((1 - x^2/a[n]^2 )*b[n]^2))
suppressWarnings(y2 <- -sqrt((1 - x^2/a[n]^2 )*b[n]^2))

#アニメーション作成
plot(x, y1, type="n",
     xlab="x", ylab="y", cex.lab=1.5, cex.axis=1.5,
     xlim=c(min(x, na.rm = T), max(x, na.rm = T)),
     ylim=c(min(x, na.rm = T), max(x, na.rm = T)))
lines(x, y1, type="l", lwd=1.2)
lines(x, y2, type="l", lwd=1.2)
plotrix::boxed.labels(x = 8, y = 8, 
                labels = paste0(" a: ", formatC(a[n], width = 2, flag = " "), "\n", 
                                " b: ", formatC(b[n], width = 2, flag = " ")),
                bg = "#b2f4f4c0", adj = 0.5, cex = 2)
plotrix::boxed.labels(x = -7, y = 8, 
                labels = expression(frac(x^2,a^2) + frac(y^2,b^2) == 1),
                bg = "#fffafa", adj = 0.5, cex = 2)

}},  movie.name = paste0("Fig_01.gif"), interval = 0.25, dpi=300,
nmax = length(a), ani.width = 600, ani.height=600, ani.type="png")

斜め楕円の作図

#斜め楕円の公式(R風)
#A*x^2 + B*x*y + C*y^2 = 1
#y = (-B*x ± sqrt((B*x)^2 -4*C*(A*x^2-1))) / 2*C 

#パラメータ設定
A1 <- c(3:-10, -10:10, 10:1, 1:10, 10:3)
B1 <- rep(3, length(A1))
C1 <- rep(3, length(A1))

B2 <- c(3:-10, -10:10, 10:1, 1:10, 10:3)
A2 <- rep(3, length(B2))
C2 <- rep(3, length(B2))

C3 <- c(3:-10, -10:10, 10:1, 1:10, 10:3)
A3 <- rep(3, length(C3))
B3 <- rep(3, length(C3))

#最終的なパラメータセット
A <- c(A1, A2, A3)
B <- c(B1, B2, B3)
C <- c(C1, C2, C3)

#アニメーション作成
animation::saveGIF({
for(n in seq_len(length(A))){
#n <- 1
#X軸
x <- seq(-15, 15, by=0.000002)

#Y軸
suppressWarnings(y1 <- (-B[n]*x + sqrt((B[n]*x)^2 -4*C[n]*(A[n]*x^2-1))) / 2*C[n])
suppressWarnings(y2 <- (-B[n]*x - sqrt((B[n]*x)^2 -4*C[n]*(A[n]*x^2-1))) / 2*C[n])

#NAを削除
x0 <- x[!is.na(y1)]
y1 <- y1[!is.na(y1)]
y2 <- y2[!is.na(y2)]

#作図
plot(x0, y1, type="n",
     xlab="x", ylab="y", cex.lab=1.5, cex.axis=1.5,
     xlim=c(min(x, na.rm = T), max(x, na.rm = T)),
     ylim=c(min(x, na.rm = T), max(x, na.rm = T)))
lines(x0, y1, type="l", lwd=0.1)
lines(x0, y2, type="l", lwd=0.1)
plotrix::boxed.labels(x = max(x, na.rm = T)*4/5, y = max(x, na.rm = T)*4/5, 
                labels = paste0(" A: ", formatC(A[n], width = 2, flag = " "), "\n", 
                                " B: ", formatC(B[n], width = 2, flag = " "), "\n",
                                " C: ", formatC(C[n], width = 2, flag = " ")),
                bg = "#b2f4f4c0", adj = 0.5, cex = 2)
plotrix::boxed.labels(x = max(x, na.rm = T)*-8/15, y = max(x, na.rm = T)*13/15, 
                labels = expression(A*x^2 + B*x*y + C*y^2 == 1),
                bg = "#fffafa", adj = 0.5, cex = 2)

}},  movie.name = paste0("Fig_02.gif"), interval = 0.2, dpi=200,
nmax = length(A), ani.width = 600, ani.height=600, ani.type="png")