Avoiding Stack overflow with Ackermann function in R

Question

I recently saw an interesting Computerphile Video about the Ackermann function and tried to recreate it in R, here's what I came up with:

Ackermann <- function(m,n){

  if (m == 0){

    return(n+1)

  } else if (m > 0 & n == 0){

    return(Ackermann(m-1,1))

  } else if (m > 0 & n > 0){

    return(Ackermann(m-1,Ackermann(m,n-1)))

  }

}

in the video, they implemented their own version of the code (in C, I think) and explained that it takes a massive amount of recursive computation for specific value pairs such as 4,1 and it took them 3 minutes to compute that value. If I try to recreate this in R with my algorithm I get a stack overflow:

Error: C stack usage  7971652 is too close to the limit

Is there a way to get the result for Ackermann(4,1) in R?

Answer 1

I think it is possible but probably quite complicated. If you write it like this (see below) it will not error out, but it will take quite some time:

sub_Ackermann1 <- function(df){
  i <- nrow(df)
  m <- df$m[i]
  n <- df$n[i]
  if (m == 0){
    r <- n+1
    df$r[i] <- r
    df_i <- df} 
  else if (m > 0 & n == 0){
    r <- NA
    m <- m-1
    n <- 1
    df_i <- df
    newrow <- data.frame(m=m,n=n,r=r)
    df_i <- rbind(df_i,newrow)} 
  else if (m > 0 & n > 0){
    r1 <- NA
    m1 <- m-1
    n1 <- NA
    df_i <- df
    newrow1 <- data.frame(m=m1,n=n1,r=r1)
    df_i <- rbind(df_i,newrow1)

    r2 <- NA
    m2 <- m
    n2 <- n-1
    newrow2 <- data.frame(m=m2,n=n2,r=r2)
    df_i <- rbind(df_i,newrow2)}

  return(df_i)
}

sub_Ackermann2 <- function(df){
  r <- df$r[nrow(df)]
  if (is.na(df$n[nrow(df)-1])){ 
    df$n[nrow(df)-1] <- r }
  else if (is.na(df$r[nrow(df)-1])){ df$r[nrow(df)-1] <- r}
  df_i <- df[-nrow(df),] 
  return(df_i)
}
Ackermann <- function(m,n){
  df <- data.frame(m=m,n=n,r=NA)
  if (m == 0){df$r <- n+1} 
  while (is.na(df$r[1])){
    if (is.na(df$r[nrow(df)])){ df <- sub_Ackermann1(df)}
    else if (is.na(df$r[1])){ df <- sub_Ackermann2(df)}
  }
  return(df$r[1])

}

It works on smaller values at least, and doesn´t crash on larger values. Maybe someone can show that this can´t work or vice versa, have ideas how to optimize it...