# Description of problem

Hi everyone,
I’ve been trying to simulate the time evolution of a neuron model while adding at each time step an arbitrary external input. This input represents the average external stimulation that my neuron population receives, therefore I decided to handle it as follows: at every time step, the value of this input represents the mean of a normal distribution, from which N values are drawn. Each neuron receives one of this values as an additive contribution to its defining differential equation, named drive.
The problem I face is that the simulation is mostly composed of nans and raises the following message:
WARNING neurongroup’s variable ‘phi’ has NaN, very large values, or encountered an error in numerical integration. This is usually a sign that an unstable or invalid integration method was chosen. [brian2.groups.group.invalid_values]
I would appreciate suggestions either on how to perform correctly the integration or on how to implement my input differently.

# Minimal code to reproduce problem

``````# configure
dtime = 100
T = 1000
N = 5000
J = 1*br.Hz
br.start_scope()
br.defaultclock.dt = dtime*br.ms

# define external drive
drive_avg = np.zeros((int(T/dt),2)) # all-zero mean to try
drive = br.TimedArray(br.stack([br.numpy.random.normal(loc=drive_avg[i,0],scale=dtime**0.5,size=N).T
for i in range(drive_avg.shape[0])],axis=1),dt=dtime*br.ms)

# define equations
eqs = '''
dphi/dt = omega_int + a_int*sin(phi) + coupl + drive(t,i)*Hz + sigma_int*xi : 1
coupl : hertz
omega_int : hertz
a_int : hertz
sigma_int : hertz**0.5
'''

# define population
neurons = br.NeuronGroup(N,eqs,threshold='sin(phi)>0.6',reset='',method='heun')
neurons.phi = 2*math.pi * br.np.random.uniform(size=N)
neurons.omega_int = 1*br.Hz
neurons.a_int = 1*br.Hz
neurons.sigma_int = 1*br.Hz**0.5
conn = br.Synapses(neurons,neurons,'coupl_post = J/N_pre*sin(phi_pre-phi_post) : hertz (summed)')
conn.connect(p=1)
activity = br.StateMonitor(neurons,'phi',record=True)

br.run(T*br.ms)
``````

# What you have aready tried

I tried reducing the variance of the normal distribution as I thought the integration error might stem from too big values.