fractal flame

ff_b.cpp

#include <Rcpp.h>
using namespace Rcpp;

// turmite function to be called from R
// [[Rcpp::export]]
NumericMatrix flame(int iter, int layers) {
  
  NumericMatrix points(iter, 3); // initially zero
  NumericMatrix coeffs(9, layers);
  
  // set coefficients
  for(int i = 0; i < 9; i++) {
    for(int j = 0; j < layers; j++) {
      coeffs(i,j) = R::runif(-1, 1);
    }
  }
  
  // initial values
  points(0, 0) = R::runif(-1, 1);
  points(0, 1) = R::runif(-1, 1);
  points(0, 2) = R::runif(-1, 1);
  
  // iterate
  int r;
  double x;
  double y;
  double z;
  double s;
  
  for(int t = 1; t < iter; t++) {
    
    r = rand() % layers; // which affine transform to use?
    
    // co-ordinates after random transform
    x = coeffs(0, r) * points(t-1, 0) + coeffs(1, r) * points(t-1, 1) + coeffs(2, r);
    y = coeffs(3, r) * points(t-1, 0) + coeffs(4, r) * points(t-1, 1) + coeffs(5, r);
    z = coeffs(6, r) * points(t-1, 0) + coeffs(7, r) * points(t-1, 1) + coeffs(8, r);
    
    // apply function to the transformed coords
    s = x*x + y*y + z*z;
    x = x/s;
    y = y/s;
    z = z/s;
    
    // store results
    points(t, 0) = x;
    points(t, 1) = y;
    points(t, 2) = (z + points(t-1, 2))/2;
    
  }
  
  return points;
}

fractal_flame_03.R

library(Rcpp)
library(ggplot2)

sourceCpp(here::here("source", "ff_b.cpp"))

# parameters
seed <- 11
iter <- 20000000
layers <- 4
bg <- "ghostwhite"
pl <- "scico::oslo"

set.seed(seed)

df <- flame(iter, layers)
df <- as.data.frame(df)
names(df) <- c("x","y","c")
df <- df[-(1:100),]

p <- ggplot(df, aes(x,y, colour = c)) + 
  geom_point(size = 1.5, alpha = .2, stroke = 0, show.legend = FALSE) + 
  theme_void() + 
  theme(plot.background = element_rect(bg, bg)) + 
  paletteer::scale_color_paletteer_c(pl)

fname <- paste0("ff_03_", seed, ".png")
fpath <- here::here("image", fname)
ggsave(fpath, p, width = 16, height = 16, dpi = 100)

library(Rcpp)
library(ggplot2)

sourceCpp(here::here("ff_b.cpp"))

# parameters
seed <- 11
iter <- 20000000
layers <- 4
bg <- "ghostwhite"
pl <- "scico::oslo"

set.seed(seed)

system.time({
df <- flame(iter, layers)
df <- df[-(1:100),]
})


# Manually scale the co-ordinates to the image size
df[,1] <- (df[,1] - min(df[,1])) / (max(df[,1]) - min(df[,1])) * 1600
df[,2] <- (df[,2] - min(df[,2])) / (max(df[,2]) - min(df[,2])) * 1600


# Manually create a vector of colours
col_idx <- as.integer((df[,3] - min(df[,3])) / (max(df[,3]) - min(df[,3])) * 255) + 1L
col <- terrain.colors(256)[col_idx]

library(cairocore)
library(cairobasic)
system.time({
  cb <- cairobasic::CairoBasic$new(width = 1600, height = 1600, antialias = FALSE)
  cb$add_circles(x=df[,1], y = df[,2], r = 1, fill = col, colour = NA)
  cb$write_png("cairo.png")
})

# About 40 seconds with cairobasic vs 4 minutes with ggsave
# I'm cheating by not aliasing.  Leaving it on doubles the time.
# still some weird behaviour with cairocore, but a quick session restart should fix things

code for rgl fwiw

rgl::par3d(windowRect = c(0, 0, 1024, 1024))
rgl::view3d(phi = 0, theta = 0, interactive = FALSE, zoom = 0.8)
system.time({
rgl::points3d(df$x, df$y, 0, color = rgb(0, 0, scales::rescale(df$c)), size = 1)
})
 user  system elapsed 
   9.40    1.65   11.53 

for completeness:

snapshot3d("rgl.png")

the file image is the dimensions of the device currently, but arbitrary render is coming soon I think

Very cool. You can also get noticeable speedups with https://github.com/exaexa/scattermore

library(Rcpp)
library(ggplot2)

sourceCpp(here::here("source", "ff_b.cpp"))

# parameters
seed <- 11
#iter <- 20000000
iter <- 5e6
layers <- 4
bg <- "ghostwhite"
pl <- "scico::oslo"

set.seed(seed)

df <- flame(iter, layers)
df <- as.data.frame(df)
names(df) <- c("x","y","c")
df <- df[-(1:100),]

system.time({
  p <- ggplot(df, aes(x,y, colour = c)) + 
    geom_point(size = 1.5, alpha = .2, stroke = 0, show.legend = FALSE) + 
    theme_void() + 
    theme(plot.background = element_rect(bg, bg)) + 
    paletteer::scale_color_paletteer_c(pl)
  
  fname <- paste0("ff_03_", seed, "_base.png")
  fpath <- here::here("image", fname)
  ggsave(fpath, p, width = 16, height = 16, dpi = 100)
})

system.time({
  p <- ggplot(df, aes(x,y, colour = c)) + 
    scattermore::geom_scattermore(pointsize = 1.5, alpha = .2, stroke = 0, show.legend = FALSE, pixels=c(1600,1600)) + 
    theme_void() + 
    theme(plot.background = element_rect(bg, bg)) + 
    paletteer::scale_color_paletteer_c(pl)
  
  fname <- paste0("ff_03_", seed, "_more.png")
  fpath <- here::here("image", fname)
  ggsave(fpath, p, width = 16, height = 16, dpi = 100)
})


dfm <- as.matrix(df[,c("x","y")])
#I am not great at color management, so there is certainly a better way... I am not including those in timing as they are sloooooow
min_col <- min(df$c)
max_col <- max(df$c)
colours <-  paletteer::paletteer_c(pl, 1000)[1 + round( (df$c - min_col) / (max_col - min_col) * 999 )]
colours_alpha <- paste0(substr(colours, 1, 7), "33") 

system.time({
  
  p <- ggplot() + scattermore::geom_scattermost(dfm, col = colours_alpha, pointsize = 1.5,  pixels=c(1600,1600))  +
    theme_void() + 
    theme(plot.background = element_rect(bg, bg)) 

  fname <- paste0("ff_03_", seed, "_most.png")
  fpath <- here::here("image", fname)
  ggsave(fpath, p, width = 16, height = 16, dpi = 100)
  
})

Thank you @coolbutuseless, @mdsumner and @martinmodrak! This is all super helpful!