Dynamic input to populate Numpy array without for loop for Monte Carlo

I've got some elegant code from other contributors for a Monte Carlo (random walk) using numpy. However, currently the 'vol,' aka standard deviation, is a supplied constant. Ideally, it should vary across steps in time. Is there a way to adapt this code to make 'vol' a function of the numpy index? I could do it with a for loop but I'm trying to keep it vectorized for performance reasons.

import numpy as np

def montecarlo_brownian(start_value, T_years, steps, vol, rate = 0, sims = 1000): 
    """ generate random walks from a starting value for Monte Carlo analysis

    Args:
        start_value (float):    starting value for random walks
        T_year (float):         number of years in question
        steps (int):            substeps within each year 
        vol (float):            annualized standard deviation (ie, implied volatility)
        rate (float):           risk free interest rate
        sims (int):             number of simulations to execute
    
    Returns:
        (np.array)  columns of random walks 
    """
    times = np.linspace(0, T_years, steps+1)
    dt = times[1] - times[0]
    
    B = np.random.normal(0, np.sqrt(dt), size=(sims, steps)).T
    S = np.exp((rate - vol ** 2/ 2) * dt + vol * B)
    S = start_value * S.cumprod(axis=0)
    return S

commenter jared requested sample of the for loop logic. I'd define a function

def get_vol(t):
  pass
  #some function that calculates vol based on how much time is left 

then use something like:

def monte_carlo_with_changing_vol(S0, r, T, num_sims, num_steps, volatility):
    dt = T / num_steps
    paths = np.zeros((num_sims, num_steps + 1))
    paths[:, 0] = S0

    for i in range(num_sims):
        for j in range(1, num_steps + 1):
            current_vol = vol(j)
            paths[i, j] = paths[i, j - 1] * np.exp((r - 0.5 * current_vol ** 2) * dt + current_vol * np.sqrt(dt) * np.random.normal(0, 1))

    return paths

If you want vol to be an array of size steps, you need to add an extra dimension for broadcasting. Because vol will have the shape (steps,) and B has the shape (steps, sims), you should add the extra dimension at the end. This can be done using None (np.newaxis will do the same thing, but that is just an alias for None), i.e. vol[:,None].

With that change, your code will look like this:

S = np.exp((rate - 0.5*vol[:,None]**2)*dt + vol[:,None]*B)