set.seed(42)

N      <- 100
years  <- 10
times  <- seq(0, years, by = 0.25)  # quarterly to keep things heavier

# Interest rates per account (0–6%/yr)
r      <- runif(N, 0.00, 0.06)

# Transfer matrix (row-stochastic-ish, mostly small)
Tmat   <- matrix(0, N, N)
idx    <- matrix(sample.int(N*N, size = N*8), ncol = 1)  # ~8 nonzeros/row on avg
Tmat[idx] <- runif(length(idx), -0.02, 0.02)             # small +/- coupling
diag(Tmat) <- diag(Tmat) - rowSums(Tmat) * 0.1           # stabilise a bit

# Seasonal deposit amplitudes and phases
A      <- runif(N, 0, 50)
phi    <- runif(N, 0, 1)

# Initial balances
B0     <- runif(N, 500, 2000)
library(odin)

net_model <- odin({
  n <- user()              # number of accounts

  r   <- user();  dim(r)   <- n
  A   <- user();  dim(A)   <- n
  phi <- user();  dim(phi) <- n
  T   <- user();  dim(T)   <- c(n, n)

  B0  <- user();  dim(B0)  <- n

  # State
  initial(B[i]) <- B0[i]
  dim(B) <- n

  # Inline the seasonal term; no intermediate D[i] needed
  deriv(B[i]) <- r[i] * B[i] + sum(j, T[i, j] * B[j]) +
                 A[i] * (1 + sin(2 * pi * (t - phi[i])))

  # Useful aggregate
  output(total) <- sum(B[i])
})

mod   <- net_model$new(n = N, r = r, A = A, phi = phi, T = Tmat, B0 = B0)

# warm run
res_o <- mod$run(times)

# benchmark (don’t include compilation)
n_runs <- 5
t_od <- system.time({
  for (k in 1:n_runs) {
    mod$run(times)
  }
})
print(t_od)

# quick visual / sanity check
plot(res_o[, "t"], res_o[, "total"], type = "l", lwd = 2,
     xlab = "Years", ylab = "Total balance",
     main = "Odin model – total balance")