impulse <- c(8e-24, 8e-24, 8e-24, 8e-24)
cume_impulse <- c(8.0e-24, 1.6e-24, 2.4e-24, 3.2e-24) #cumulative total
time_steps <- c(1.132, 6.136, 11.140, 16.144) #time step at which the increase occurs
bar_data <- data.frame(time_steps, impulse, cume_impulse)

library(ggplot2)
library(gganimate)
library(dplyr)
library(gifski)
library(tidyverse)
library(magick)
plots <- bar_data %>% 
  ggplot(aes(x = x_pos, y = cume_impulse)) + geom_col()
#label the animation for the bar plot as bar_anim
bar_anim <- plots + transition_time(time_steps)
#set the parameters and display the bar plot animation
bar_gif <- animate(bar_anim, height = 500, width = 800, fps = 20, duration = 20)
#save the bar plot to a gif file
anim_save("bar_impulse.gif", bar_gif)
#the following two lines are required before the image_append function will process
mgif_particle <- image_read(particle_gif) 
mgif_bar <- image_read(bar_gif)
#join the first fraomes of the gifs
paired_gif <- image_append(c(mgif_particle[1], mgif_bar[1]))
for(i in 2:317){
  combined <- image_append(c(mgif_particle[i], mgif_bar[i]))
  paired_gif <- c(paired_gif, combined)
}
anim_save("bar_plus_particle.gif", paired_gif)

m <- 1e-26    #Mass of a particle in kg
vx_i <- 400   #initial x velocity in m/s
vy_i <- 200   #initial y velocity in m/s
#creating variables for the initial coordinates of the particle
x_i <- round(runif(n = 1, min = 0,max = 1), 3)
y_i <- round(runif(n = 1, min = 0,max = 1), 3)
#initialise vectors that will be added to a data frame
x_pos <- round(seq(x_i, (x_i + 17), 0.001), 3) #17 was an arbitrary choice to allow for a few bounces on the right wall
y_pos <- round(seq(y_i, (y_i + 17), 0.001), 3) 
time_steps <- seq(0, 17.001, 0.001)
#loop to detect x direction wall contact
xv = 0.0004
for (i in 1 + seq_along(x_pos)) {
  if (x_pos[i-1] > 1 | x_pos[i-1] < 0 ){xv = xv*(-1)}
  x_pos[[i]] <- round(x_pos[[i-1]] + xv, 4)
  }
#loop for y direction wall contact
yv = 0.0002
for (i in 1 + seq_along(y_pos)) {
  if (y_pos[i-1] > 1 | y_pos[i-1] < 0 ){yv = yv*(-1)}
  y_pos[[i]] <- round(y_pos[[i-1]] + yv, 4)
  }
#creating a data frame with the particle's kinematic information (this line must come after the loops)
pos_data <- data.frame(x_pos, y_pos, time_steps)
graph1 = pos_data %>% ggplot(aes(x_pos, y_pos)) + 
  geom_point(colour = "red", size = 2) + 
  xlab(NULL) + ylab(NULL) + xlim(0, 1) + ylim(0, 1)
#initialising vector for x-momentum. Same length as x_pos with arbitrary values.
px <- rep(xv,length(x_pos)) 
#add a vector to pos_data called x_dif to pos_data that represents the change in x position between successive time steps
pos_data$x_dif <- c(0, diff(pos_data$x_pos))
# add a vector to pos_data called x_turn which will help detect changes in direction
pos_data$x_turn <- c(0, diff(sign(pos_data$x_dif)))
# add a cumulative total column for impulse on right wall called cume_impulse
pos_data$cume_impulse <- cumsum(pos_data$impulse)
#add a vector to pos_data which will include the impulse magnitude on the right wall at each time step
pos_data <- pos_data %>% mutate(impulse = case_when(abs(pos_data$x_turn) == 2 ~ 2*(vx_i*m),
                           abs(pos_data$x_turn) != 2 ~ 0))
graph1_animation = graph1 + transition_time(time_steps) +
  labs(subtitle = "{frame_time}") + ease_aes('linear') +
  shadow_wake(wake_length = 0.05) #animation code
particle_gif <- animate(graph1_animation, height = 500, width = 800, nframes = 317)
anim_save("particle_bounce.gif", particle_gif) #saves animation as a gif file
#Filter the data to only include frames that include impulse, and x_pos = 1 (right wall)
bar_data <- filter(pos_data,impulse!=0 & x_pos == 1)
#new column in bar_data called "cume_impulse" - cumulative total of the impulse on the right wall
bar_data["cume_impulse"] <- cumsum(bar_data$impulse)

## initial position for particle
x_i <- round(runif(n = 1, min = 0,max = 1), 3)
y_i <- round(runif(n = 1, min = 0,max = 1), 3)
## more initial conditions
m <- 1e-26    #Mass of a particle in kg
vx_i <- 400   #initial x velocity in m/s
vy_i <- 200   #initial y velocity in m/s
#initialise vectors that will be added to a data frame
x_pos <- round(seq(x_i, (x_i + 17), 0.001), 3)
y_pos <- round(seq(y_i, (y_i + 17), 0.001), 3) 
time_steps <- seq(0, 17.001, 0.001)
#xv_column <- seq(0, 17.001, 0.001) # don't think this is necessary
#loop for x direction wall contact
xv = 0.0004
for (i in 1 + seq_along(x_pos)) {
  if (x_pos[i-1] > 1 | x_pos[i-1] < 0 ){xv = xv*(-1)}
  x_pos[[i]] <- round(x_pos[[i-1]] + xv, 4)
}
#loop for y direction wall contact
yv = 0.0002
for (i in 1 + seq_along(y_pos)) {
  if (y_pos[i-1] > 1 | y_pos[i-1] < 0 ){yv = yv*(-1)}
  y_pos[[i]] <- round(y_pos[[i-1]] + yv, 4)
}
#creating a data frame with the particle's kinematic information (this line must come after the loops)
pos_data <- data.frame(x_pos, y_pos, time_steps)
#initialising vector for x-momentum. Same length as x_pos with arbitrary values.
px <- rep(xv,length(x_pos)) 
#add a vector to pos_data called x_dif to pos_data that represents the change in x position between successive time steps
pos_data$x_dif <- c(0, diff(pos_data$x_pos))
# add a vector to pos_data called x_turn which will help detect changes in direction
pos_data$x_turn <- c(0, diff(sign(pos_data$x_dif)))
#add a vector to pos_data which will include the impulse magnitude on the right wall at each time step
pos_data <- pos_data %>% mutate(impulse = case_when(abs(pos_data$x_turn) == 2 ~ 2*(vx_i*m),
                                                    abs(pos_data$x_turn) != 2 ~ 0))
# add a cumulative total column for impulse on right wall called cume_impulse
pos_data$cume_impulse <- cumsum(pos_data$impulse)
## ------------------------------------ tjebo's ------------------------------------------------------
## --- I may have misunderstood which data you used.Did you use the the length 4 vectors below for the bar_data frame? 
## ---- I presumed not, so I made the bar_data frame similar length to pos_data.
## using your pos_data frame and impulse vector
## I have changed the order of cume_impulse, because it somehow made more sense
#cume_impulse <- c(1.6e-24, 2.4e-24, 3.2e-24, 8.0e-24)
## time steps need to be character, for the later merge
#time_steps <- as.character(c(1.132, 6.136, 11.140, 16.144))
bar_data <- data.frame(pos_data$time_steps, pos_data$impulse, pos_data$cume_impulse)
df1 <- pos_data %>% mutate(time_steps = as.character(time_steps))
# --- data frames wouldn't merge because of the automatic name change to the columns. This was how I fixed that. Unsure why this didn't effect you
names(bar_data)[names(bar_data)=="pos_data.time_steps"] <- "time_steps"     
names(bar_data)[names(bar_data)=="pos_data.cume_impulse"] <- "cume_impulse"
names(bar_data)[names(bar_data)=="pos_data.impulse"] <- "impulse"
# --- I fragmented the code to help me understand the individual commands, and also to locate where it was breaking down
df2 <- pos_data %>% left_join(bar_data, by = "time_steps")
df3 <- pos_data %>% bind_rows(., .)
df4 <- df3 %>% mutate(
  time_steps = as.numeric(time_steps), 
  panel = rep(c("bar", "dot"), each = nrow(pos_data)), 
  ## removing x and y for the bar panel
  across(c(x_pos, y_pos), ~ifelse(panel == "bar", NA, .x))
)
## --- I wasn't quite sure what this command was for
df5 <- df4 %>% fill(cume_impulse)
## --- I presumed this converted the NA's into zeros in cume_impulse?
df5$cume_impulse[is.na(df5$cume_impulse)] <- 0
## --- Looks like this converts the cume_impulse values that correspond to dot as NA. Unsure why it's necessary though.
df5$cume_impulse[df5$panel == "dot"] <- NA
set.seed(1)
times <- sort(sample(unique(df5$time_steps), 100))
df_frac <- df5[df5$time_steps %in% times, ]
## dummy data for facet range
bar_range <- data.frame(x= .5, cume_impulse = range(df5$cume_impulse, na.rm = T), panel = "bar")
## separate out the aesthetics to the respective layers
p <-
  ggplot(df_frac) +
  geom_blank(data = bar_range, aes(x = .5, y = cume_impulse)) +
  ## need position = "identity"
  geom_col(aes(x = .5, y = cume_impulse), na.rm = T, position = "identity") +
  geom_point(aes(x_pos, y_pos), na.rm = T) +
  facet_wrap(~panel, scales = "free_y") +
  scale_x_continuous(expand = c(0,0))
p_anim <- p + transition_time(time_steps)
# this is not random = there are 100 unique time steps, and now 100 frames 
animate(p_anim, fps = 20, duration = 5)

suppressMessages({
  library(ggplot2)
  library(dplyr)
  library(tidyr)
  library(gganimate)
})
## using your pos_data frame and impulse vector
## I have changed the order of cume_impulse, because it somehow made more sense
cume_impulse <- c(1.6e-24, 2.4e-24, 3.2e-24, 8.0e-24)
## time steps need to be character, for the later merge
time_steps <- as.character(c(1.132, 6.136, 11.140, 16.144))
bar_data <- data.frame(time_steps, impulse, cume_impulse)
df <-
  pos_data %>%
  ## you need time steps as character for a correct merging - probably a floating point issue
  mutate(time_steps = as.character(time_steps)) %>%
  left_join(bar_data, by = "time_steps")  %>%
  bind_rows(., .) %>%
  ## adding a facetting variable
  mutate(
    time_steps = as.numeric(time_steps), 
    panel = rep(c("bar", "dot"), each = nrow(pos_data)), 
  ## removing x and y for the bar panel
    across(c(x_pos, y_pos), ~ifelse(panel == "bar", NA, .x))
  ) %>%
## giving values to all time steps
  fill(cume_impulse) 
df$cume_impulse[is.na(df$cume_impulse)] <- 0
## removing y for dot panel 
df$cume_impulse[df$panel == "dot"] <- NA
## I am reducing the data to a mere fraction, because you get the same plot with way
## less computation (also, otherwise my computer is crashing). 
## I do this after merging etc so not to accidentally remove the bar values
## I'm deliberately selecting 100 unique time steps for the animation (so that frame number is equivalent to time stamps)
set.seed(1)
times <- sort(sample(unique(df$time_steps), 100))
df_frac <- df[df$time_steps %in% times, ]
## dummy data for facet range
bar_range <- data.frame(x= .5, cume_impulse = range(df$cume_impulse, na.rm = T), panel = "bar")
## separate out the aesthetics to the respective layers
p <-
  ggplot(df_frac) +
  geom_blank(data = bar_range, aes(x = .5, y = cume_impulse)) +
  ## need position = "identity"
  geom_col(aes(x = .5, y = cume_impulse), na.rm = T, position = "identity") +
  geom_point(aes(x_pos, y_pos), na.rm = T) +
  facet_wrap(~panel, scales = "free_y") +
  scale_x_continuous(expand = c(0,0))
p_anim <- p + transition_time(time_steps)
## this is not random = there are 100 unique time steps, and now 100 frames 
animate(p_anim, fps = 20, duration = 5)