一个项目的一部分是模拟气体粒子对墙壁的冲击，测量随着时间推移的作用力。到目前为止，我已经成功地完成了一个粒子的模拟，在考虑如何添加更多粒子时，我感觉到有一种更有效的方法来构建代码。
我的想法基本上是对第一个粒子做同样的事情（用相关数据建立一个数据框），然后为每个要添加到图中的新粒子添加另一个geom_point。
我会很感激一些建议，如何添加更多的粒子到这个模拟，而不必有一个脚本，其中包含了很多类似的段落。
以下是到目前为止我对单个粒子的理解：

library(tidyverse)
library(magick)
library(gganimate)
library(stringi)
library(ggtext)

particle_mass <- 4.65 #(x 10^-26) 
## A constant that dictates the number of decimal places that the data will be rounded to. Having this set as a constant avoids having to manually adjust a number which appears multiples times throughout the code if you wanted to change the number of decimal places for any reason.
rnd <- 4 
## initial particle coordinates 
x_i <- round(runif(n = 1, min = 0,max = 1), rnd) 
y_i <- round(runif(n = 1, min = 0,max = 1), rnd)

x_velocity <- 400 ## A constant for the x velocity in m/s
y_velocity <- 200 ##A  constant for the y velocity in m/s
time_step <- 0.001 ## A constant for the duration of a single time step in seconds
frame_total <- 200 ## A constant for the number of frames that the animation will be split into (two additional frames are added on at the end)
rate_factor <- 35 ## A factor that controls the playback rate of the animation. Increase this number and the animation will be slowed 

## The x and y displacement components per time step with a factor to slow down the animation to a more easily perceptible speed.
x_jump <- (x_velocity/rate_factor)/frame_total 
y_jump <- (y_velocity/rate_factor)/frame_total 

## Vectors of equal length for x position, y position, and corresponding time steps are initialized.
x_pos <- round(seq(x_i, (x_i + (frame_total*x_jump))-x_jump, x_jump), rnd)
y_pos <- round(seq(y_i, (y_i + (frame_total*y_jump))-y_jump, y_jump), rnd) 
times <- seq(0, (frame_total*time_step)-time_step, time_step)

## A loop that adjusts the direction of the x velocity vector when x = 1 or x = 0 to simulate wall contact
for (i in 2:length(x_pos)) {
  if (x_pos[i-1] > 1 | x_pos[i-1] < 0 ){x_jump = x_jump*(-1)}
  x_pos[[i]] <- round(x_pos[[i-1]] + x_jump, rnd)
}

##  A loop that adjusts the direction of the y velocity vector when y = 1 or y = 0 to simulate wall contact
for (i in 2:length(y_pos)) {
  if (y_pos[i-1] > 1 | y_pos[i-1] < 0 ){y_jump = y_jump*(-1)}
  y_pos[[i]] <- round(y_pos[[i-1]] + y_jump, rnd)
}

bounce <- data.frame(times, x_pos, y_pos) ##time and coordinate position information
bounce$x_dif <- c(0, diff(bounce$x_pos))## change in x position between time steps
bounce$x_turn <- c(0, diff(sign(bounce$x_dif))) ## helps detect changes in direction

## impulse magnitude experienced by the right wall at each time step
bounce <- bounce %>% mutate(impulse_right = case_when((bounce$x_turn) == -2 ~ 2*(x_velocity*particle_mass),
                                                      (bounce$x_turn) != -2 ~ 0))

bounce$cume_impulse_right <- cumsum(bounce$impulse_right) #cumulative impulse on the right wall 
## static plots of the moving particle
bounce_plot <- ggplot(bounce) +
  geom_point(aes(x_pos, y_pos), color = "azure4", size = 1) +
  labs(title = "Single gas particle model") + 
  ylab("") + xlab("") +
  theme_classic() +
  ## axis lines and tick marks
  theme(axis.line = element_blank(),
        axis.ticks.x = element_blank(),
        axis.text.x = element_blank(), axis.ticks.y = element_blank(),
        axis.text.y = element_blank()) +
  ## plot title aesthetics
  theme(plot.title = element_text(family = "Arial", face = "bold", hjust = 0.5, 
                                  size = 9, color = "gray38")) +
  ## The vertical and horizontal lines create a border around the 1 x 1 box that contains the particle
  geom_vline(xintercept=0, color="snow3", size=2) + 
  geom_hline(yintercept=0, color="snow3", size=2) + 
  geom_hline(yintercept=1, color="snow3", size=2) +
  geom_vline(xintercept=1, color="#CC6666", size=1) 

bounce_anim <- bounce_plot + transition_time(times)

## The desired number of frames and duration of the animation are set within the 'animate' function.
bounce_gif <- animate(bounce_anim, width = 3, height = 3, units = "in", res = 200,
                      nframes = frame_total, renderer = magick_renderer())

我的计划是为每个新粒子创建一个新的数据框，然后使用新数据在ggplot函数下添加一条新的geom_point线。我担心如果我添加5个，10个或更多粒子，代码会很长。有没有一种方法可以在模拟中添加多个粒子而无需编写这么多代码？任何帮助都将不胜感激。