A more streamlined approach than nested for loops to fill an array in R

I'm trying to write a function in R that models population sizes (of different age brackets) across a network of connected populations. As a simplified example, imagine I want to predict how many juveniles and adults there are in 2 separate populations, over 10 discrete time steps:

age <- c("juvenile","adult")
pop<-c(1:2)
t <- c(1:10)

I make an empty 3D array to fill, and set starting (t=1) population sizes of 2 and 4.

a1 <- array(0, dim=c(10,2,2), dimnames = list(t, age,pop))
a1[1,2,1:2] <- c(2, 4)

And then I've written a function to fill this array from t=2 to t=10:

f1 <- function(a1){
  for(t in 2:10){
    for(p in 1:2){
    a1[t,1,1] <- a1[t-1,2,1] * 2
    a1[t,2,1] <- a1[t-1,1,1] + (a1[t-1,2,1] * 0.5)
    
    a1[t,1,2] <- a1[t-1,2,2] * 3
    a1[t,2,2] <- a1[t-1,1,2] + (a1[t-1,2,2] * 0.25)
    }
  }
  return(data.frame(a1))
}

f1(a1)

Between each time step:

1. juveniles become adults,
2. adults produce more juveniles (either 2 or 3 per adult, depending on population), and
3. a proportion of adults (either 0.5 or 0.25, depending on pop) survive to the next time step.

I quite like this methodical approach, but in reality my array has many more dimensions, all of different lengths, that are all inter-connected. In other words, a whole lot of nested for loops and excessive array co-ordinates (it gets tricky to remember what something like a1[t,3,1,2,1] is actually referring to).

If anyone has any suggestions on where to start for a more streamlined approach, that would be greatly appreciated. I keep hearing how apply functions are more concise than for loops, but I've only been able to find examples on how to use apply to summarise arrays, rather than use it to fill them.