/*
* main.cpp
*
* This file is part of the refactored Izhikevich polychronization model application.
*
* This source code is based on the poly_spnet.cpp and spnet.cpp source code available at
* https://www.izhikevich.org/publications/spnet.htm.
*
* Copyright (C) 2018, Author: G. Trensch
*
* The refactored Izhikevich polychronization model application is free software:
* you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 2 of the License, or
* (at your option) any later version.
*
* It is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this application. If not, see .
*
*/
#include "params.h"
#if( __RUN_REFACTORED_VERSION__ )
#include
#include
#include
#include
#include "globals.h"
#include "utils.h"
// = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =
// = C H E C K P A R A M E T E R S E T T I N G S
// = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =
#ifndef SIM_TIME
#error SIM_TIME not defined!
#endif
#ifdef GENERATE_NETWORK_FROM_EXTERNAL_DATA
#if( GENERATE_NETWORK_FROM_EXTERNAL_DATA )
#ifndef INFILE_CONNECTION_MATRIX
#error INFILE_CONNECTION_MATRIX not defined!
#endif
#ifndef INFILE_DELAY_MATRIX
#error INFILE_DELAY_MATRIX not defined!
#endif
#ifndef INFILE_WEIGHT_MATRIX
#error INFILE_WEIGHT_MATRIX not defined!
#endif
#else
#ifndef OUTFILE_CONNECTIONS
#error OUTFILE_CONNECTIONS not defined!
#endif
#ifndef OUTFILE_DELAYS
#error OUTFILE_DELAYS not defined!
#endif
#ifndef OUTFILE_STIMULUS
#error OUTFILE_STIMULUS not defined!
#endif
#endif
#else
#error GENERATE_NETWORK_FROM_EXTERNAL_DATA not defined!
#endif
#ifdef USE_EXTERNAL_STIMULUS
#if( USE_EXTERNAL_STIMULUS )
#ifndef INFILE_STIMULUS
#error INFILE_STIMULUS not defined!
#endif
#endif
#else
#error USE_EXTERNAL_STIMULUS not defined!
#endif
#ifndef USE_STDP
#error USE_STDP not defined!
#endif
#ifndef OUTFILE_FIRINGS
#error OUTFILE_FIRINGS not defined!
#endif
// = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =
// = F O R W A R D D E C L A R A T I O N S
// = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =
void InitializeNetwork();
void InitializeSimulation();
void FinalizeSimulation();
// = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =
// = M A I N E N T R Y
// = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =
int main() {
printf( "\n\n%s\n\n", PARAM_INFO_STRING );
printf( "[INFO] Running simulation for %d seconds\n", SIM_TIME );
#if( ODE_SOLVER_REFINEMENT )
for( int i = 0; i < NUM_TOTAL_NEURONS; ++i ) {
spike[i] = 0;
}
#endif
InitializeSimulation();
InitializeNetwork();
for( int simTimeSecond = 0; simTimeSecond < SIM_TIME; ++simTimeSecond ) { // simulation loop [seconds]
for( int t = 0; t < 1000; ++t ) { // simulation loop [milliseconds]
// clear I_ext[] and select a random neuron for input
for( int i = 0; i < NUM_TOTAL_NEURONS; ++i ) {
I_ext[i] = 0.0;
}
// REFACTOR COMMENT: the original implementation does not allow networks sizes < 1000 neurons
// for( int idx = 0; idx < NUM_TOTAL_NEURONS / 1000; ++idx ) { ... }
#if( USE_EXTERNAL_STIMULUS )
int inputNeuron = GetNextExternalStimulusFromFile(INFILE_STIMULUS, simTimeSecond, t);
#else
int inputNeuron = GET_RANDOM_INT( NUM_TOTAL_NEURONS );
#endif
if( inputNeuron > 0 && inputNeuron != RC_NOID) {
I_ext[inputNeuron] = I_EXT;
}
#ifdef OUTFILE_STIMULUS
RecordRandomStimulusToFile( OUTFILE_STIMULUS, simTimeSecond, t, inputNeuron );
#endif
#if( ODE_SOLVER_REFINEMENT )
for( int n = 0; n < NUM_TOTAL_NEURONS; ++n ) {
if( spike[n] > 0 ) { // Does neuron n has fired?
if( spike[n] > 1 ) printf( "[INFO] More than 1 spike in simulation interval.\n" );
spike[n] = 0;
LTP[n][t + MAX_SYNAPSE_DELAY] = 0.1;
LTD[n] = 0.12;
for( int idx = 0; idx < numPreSynapticNeuronsOfTarget[n]; ++idx ) { // LTP: pre-synaptic spike precedes post-synaptic spike
int preSynNeuron = listOfPresynapticNeurons[n][idx];
int preSynNeuron_correspondingDelay = listOfPresynapticDelays[n][idx];
*listOfPointersToSynapticWeights_derivatives[n][idx] += LTP[preSynNeuron][t + MAX_SYNAPSE_DELAY - preSynNeuron_correspondingDelay - 1];
}
firings[numFirings][TIME] = t;
firings[numFirings][NEURON] = n;
numFirings++;
if( numFirings == MAX_NUM_FIRINGS) {
printf( "[WARNING] Two many spikes at t = %d (ignoring all)\n", t );
numFirings = 1;
}
}
}
#else
for( int n= 0; n < NUM_TOTAL_NEURONS; ++n ) {
if( v[n] >= RS_FS_THR ) { // Does neuron n has fired?
v[n] = RS_FS_C; // threshold dynamics
u[n] += d[n];
LTP[n][t + MAX_SYNAPSE_DELAY] = 0.1;
LTD[n] = 0.12;
for ( int idx = 0; idx < numPreSynapticNeuronsOfTarget[n]; ++idx ) { // LTP: pre-synaptic spike precedes post-synaptic spike
int preSynNeuron = listOfPresynapticNeurons[n][idx];
int preSynNeuron_correspondingDelay = listOfPresynapticDelays[n][idx];
*listOfPointersToSynapticWeights_derivatives[n][idx] += LTP[preSynNeuron][t + MAX_SYNAPSE_DELAY - preSynNeuron_correspondingDelay - 1];
}
firings[numFirings][TIME] = t;
firings[numFirings][NEURON] = n;
numFirings++;
if( numFirings == MAX_NUM_FIRINGS ) {
printf( "[INFO] Two many spikes at t = %d (ignoring all).\n", t );
numFirings = 1;
}
}
}
#endif
// Go back through firings as far as MAX_SYNAPSE_DELAY lasts. Calculate I_ext for the next time step.
// Take the delays of the synapses into account. I_ext applies in the next simulation step.
//
// | | |
// | | | | | | | | <-- spikes
// ---+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--> simulation time
// -5 -4 -3 -2 -1 0 1 2 3 4 5 6 7 8 9 10 11 12
// t
//
// | |
// +--------------+ MAX_SYNAPSE_DELAY (e.g., 5)
int idx = numFirings - 1; // array index starts with 0
while( t - firings[idx][TIME] < MAX_SYNAPSE_DELAY ) {
int neuronFired = firings[idx][NEURON];
int timeSinceNeuronFired = t - firings[idx][TIME];
int numEntriesOfDelay = numEntriesPerDelay[neuronFired][timeSinceNeuronFired];
for( int i = 0; i < numEntriesOfDelay; ++i ) {
int s = listOfSynsByNeuronAndDelay[neuronFired][timeSinceNeuronFired][i];
int postSynNeuron = matrixOfPostSynapticNeurons[neuronFired][s];
I_ext[postSynNeuron] += matrixOfSynapticWeights[neuronFired][s];
if( neuronFired < NUM_EXCITATORY_NEURONS) { // LTD: this spike is before postsynaptic spikes
matrixOfSynapticWeights_derivatives[neuronFired][s] -= LTD[postSynNeuron];
}
}
idx--;
}
// = = = = = = = = = = = = = = = = = = = = = =
// = UPDATE NETWORK STATE
// = = = = = = = = = = = = = = = = = = = = = =
for( int n = 0; n < NUM_TOTAL_NEURONS; ++n ) {
#if( ODE_SOLVER_REFINEMENT )
for( int intCycles = 0; intCycles < ODE_SOLVER_STEPS; ++intCycles ) {
v[n] += (1.0 / ODE_SOLVER_STEPS) * ((0.04 * v[n] + 5) * v[n] + 140 - u[n] + I_ext[n]);
u[n] += (1.0 / ODE_SOLVER_STEPS) * (a[n] * (RS_FS_B * v[n] - u[n]));
#if( LOG_MIN_MAX_V_U )
// log min and max values of v(t) and u(t)
if( v[n] > v_max ) v_max = v[n];
if( v[n] < v_min ) v_min = v[n];
if( u[n] > u_max ) u_max = u[n];
if( u[n] < u_min ) u_min = u[n];
#endif
// threshold detection for exact integration
if( v[n] >= 30.0 ) {
v[n] = RS_FS_C;
u[n] += d[n];
spike[n]++; // remember the spike event, which is aligned to next grid point
}
}
#else // original implementation
v[n] += 0.5 * ((0.04 * v[n] + 5) * v[n] + 140 - u[n] + I_ext[n]); // for numerical stability
v[n] += 0.5 * ((0.04 * v[n] + 5) * v[n] + 140 - u[n] + I_ext[n]); // time step is 0.5 ms
u[n] += a[n] * (RS_FS_B * v[n] - u[n]);
#if( LOG_MIN_MAX_V_U ) // log min and max values of v(t) and u(t)
if(v[n] > v_max) v_max = v[n];
if(v[n] < v_min) v_min = v[n];
if(u[n] > u_max) u_max = u[n];
if(u[n] < u_min) u_min = u[n];
#endif
#endif
LTP[n][t + MAX_SYNAPSE_DELAY + 1] = 0.95 * LTP[n][t + MAX_SYNAPSE_DELAY];
LTD[n] *= 0.95;
}
} // end of simulation loop [milliseconds]
// every minute: report on the average firing rate of the excitatory and inhibitory population
if( simTimeSecond > 0 && (simTimeSecond + 1) % 60 == 0 ) {
int spikesTotalExcitatory = 0;
int spikesTotalInhibitory = 0;
for( int idx = 0; idx < numFirings; ++idx ) {
if( firings[idx][NEURON] < NUM_EXCITATORY_NEURONS) {
spikesTotalExcitatory++;
}
else {
spikesTotalInhibitory++;
}
}
printf( "[INFO] Simulation time: %d seconds\n", simTimeSecond + 1 );
printf( "[INFO] ... Firing rate (excitatory) = %f\n",
(double) spikesTotalExcitatory / (double) NUM_EXCITATORY_NEURONS);
printf( "[INFO] ... Firing rate (inhibitory) = %f\n",
(double) spikesTotalInhibitory / (double) NUM_INHIBITORY_NEURONS);
#if(LOG_MIN_MAX_V_U)
printf( "[INFO] ... v_max = %f\n", v_max );
printf( "[INFO] ... v_min = %f\n", v_min );
printf( "[INFO] ... u_max = %f\n", u_max );
printf( "[INFO] ... u_min = %f\n", u_min );
#endif
}
// = = = = = = = = = = = = = = = = = = = = = =
// = SAVE FIVE SELECTED NETWORK STATES
// = = = = = = = = = = = = = = = = = = = = = =
#ifdef SELECTED_STATE_1_AFTER_N_SECONDS
#ifndef OUTFILE_STATE_1_WEIGHT_MATRIX
#error OUTFILE_STATE_1_WEIGHT_MATRIX not defined!
#endif
if( simTimeSecond == SELECTED_STATE_1_AFTER_N_SECONDS - 1 ) {
ExportWeightMatrixToFile( OUTFILE_STATE_1_WEIGHT_MATRIX );
}
#endif
#ifdef SELECTED_STATE_2_AFTER_N_SECONDS
#ifndef OUTFILE_STATE_2_WEIGHT_MATRIX
#error OUTFILE_STATE_2_WEIGHT_MATRIX not defined!
#endif
if( simTimeSecond == SELECTED_STATE_2_AFTER_N_SECONDS - 1 ) {
ExportWeightMatrixToFile( OUTFILE_STATE_2_WEIGHT_MATRIX );
}
#endif
#ifdef SELECTED_STATE_3_AFTER_N_SECONDS
#ifndef OUTFILE_STATE_3_WEIGHT_MATRIX
#error OUTFILE_STATE_3_WEIGHT_MATRIX not defined!
#endif
if( simTimeSecond == SELECTED_STATE_3_AFTER_N_SECONDS - 1 ) {
ExportWeightMatrixToFile( OUTFILE_STATE_3_WEIGHT_MATRIX );
}
#endif
#ifdef SELECTED_STATE_4_AFTER_N_SECONDS
#ifndef OUTFILE_STATE_4_WEIGHT_MATRIX
#error OUTFILE_STATE_4_WEIGHT_MATRIX not defined!
#endif
if( simTimeSecond == SELECTED_STATE_4_AFTER_N_SECONDS - 1 ) {
ExportWeightMatrixToFile( OUTFILE_STATE_4_WEIGHT_MATRIX );
}
#endif
#ifdef SELECTED_STATE_5_AFTER_N_SECONDS
#ifndef OUTFILE_STATE_5_WEIGHT_MATRIX
#error OUTFILE_STATE_5_WEIGHT_MATRIX not defined!
#endif
if( simTimeSecond == SELECTED_STATE_5_AFTER_N_SECONDS - 1 ) {
ExportWeightMatrixToFile( OUTFILE_STATE_5_WEIGHT_MATRIX );
}
#endif
RecordNetworkActivityToFile( OUTFILE_FIRINGS, simTimeSecond, numFirings );
// = = = = = = = = = = = = = = = = = = = = = =
// = PREPARE NEXT SIMULATION TIME STEP
// = = = = = = = = = = = = = = = = = = = = = =
for( int n = 0; n < NUM_TOTAL_NEURONS; ++n ) {
for( int i = 0; i < MAX_SYNAPSE_DELAY + 1; ++i ) {
LTP[n][i] = LTP[n][1000 + i];
}
}
int k = numFirings - 1;
while( 1000 - firings[k][TIME] < MAX_SYNAPSE_DELAY) {
k--;
}
for( int i = 1; i < numFirings - k; ++i ) {
firings[i][TIME] = firings[k + i][TIME] - 1000;
firings[i][NEURON] = firings[k + i][NEURON];
}
numFirings = numFirings - k;
#if( USE_STDP )
// = = = = = = = = = = = = = = = = = = = = = =
// = UPDATE SYNAPTIC WEIGHTS
// = = = = = = = = = = = = = = = = = = = = = =
// only excitatory connections are modified
for( int n = 0; n < NUM_EXCITATORY_NEURONS; ++n ) {
for( int s = 0; s < NUM_SYNAPSES_PER_NEURON; ++s ) {
// REFACTOR COMMENT:
// The following two lines have a different order in spnet.cpp.
matrixOfSynapticWeights_derivatives[n][s] *= 0.9;
matrixOfSynapticWeights[n][s] += 0.01 + matrixOfSynapticWeights_derivatives[n][s];
if( matrixOfSynapticWeights[n][s] > MAX_SYNAPTIC_STRENGTH) {
matrixOfSynapticWeights[n][s] = MAX_SYNAPTIC_STRENGTH;
}
if( matrixOfSynapticWeights[n][s] < 0 ) {
matrixOfSynapticWeights[n][s] = 0.0;
}
}
}
#endif
} // end of simulation loop [seconds]
FinalizeSimulation();
printf( "\n\nSimulation terminated normally! \n\n" );
return (0);
}
// = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =
// = Initialize the polychronization network
// = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =
void InitializeNetwork() {
// initialize izhikevich neuron parameters
for( int n = 0; n < NUM_TOTAL_NEURONS; ++n ) {
if( n < NUM_EXCITATORY_NEURONS ) { // excitatory neurons, RS type
a[n] = RS_A;
d[n] = RS_D;
v[n] = RS_V_INIT;
u[n] = RS_U_INIT;
}
else { // inhibitory neurons, FS type
a[n] = FS_A;
d[n] = FS_D;
v[n] = FS_V_INIT;
u[n] = FS_U_INIT;
}
}
// = = = = = = = = = = = = = = = = = = = = = =
// = CONNECTION MATRIX
// = = = = = = = = = = = = = = = = = = = = = =
#if( GENERATE_NETWORK_FROM_EXTERNAL_DATA )
ImportConnectionMatrixFromFile( INFILE_CONNECTION_MATRIX );
#else
// fill matrixOfPostSynapticNeurons[n][n] with random target neurons
// avoid self-assignments and multiple connections
for( int n = 0; n < NUM_TOTAL_NEURONS; ++n ) {
for( int s = 0; s < NUM_SYNAPSES_PER_NEURON; ++s ) {
bool duplicateTarget = false;
bool selfAssigned = false;
int randomTargetNeuron = 0;
do {
duplicateTarget = false;
selfAssigned = false;
if( n < NUM_EXCITATORY_NEURONS ) { // excitatory neurons can connect to all neurons
randomTargetNeuron = GET_RANDOM_INT( NUM_TOTAL_NEURONS );
}
else { // inhibitory neurons can connect to excitatory neurons only
randomTargetNeuron = GET_RANDOM_INT( NUM_EXCITATORY_NEURONS );
}
if( randomTargetNeuron == n ) {
selfAssigned = true;
}
for( int i = 0; i < s; ++i ) {
if( matrixOfPostSynapticNeurons[n][i] == randomTargetNeuron ) {
duplicateTarget = true;
}
}
} while( duplicateTarget || selfAssigned );
matrixOfPostSynapticNeurons[n][s] = randomTargetNeuron;
}
}
#endif
#ifdef OUTFILE_CONNECTIONS
ExportConnectionMatrixToFile( OUTFILE_CONNECTIONS );
#endif
#if __DEBUG__
PrintMatrixOfPostSynapticNeurons();
#endif
// = = = = = = = = = = = = = = = = = = = = = =
// = WEIGHT MATRIX
// = = = = = = = = = = = = = = = = = = = = = =
#if(GENERATE_NETWORK_FROM_EXTERNAL_DATA)
ImportWeightMatrixFromFile(INFILE_WEIGHT_MATRIX);
#else
for( int n = 0; n < NUM_TOTAL_NEURONS; ++n ) {
for( int s = 0; s < NUM_SYNAPSES_PER_NEURON; ++s ) {
if( n < NUM_EXCITATORY_NEURONS) {
matrixOfSynapticWeights[n][s] = INIT_EXC_SYNAPTIC_WEIGHT;
}
else {
matrixOfSynapticWeights[n][s] = INIT_INH_SYNAPTIC_WEIGHT;
}
matrixOfSynapticWeights_derivatives[n][s] = 0.0;
}
}
#endif
#ifdef OUTFILE_WEIGHTS_INITIAL
ExportWeightMatrixToFile(OUTFILE_WEIGHTS_INITIAL);
#endif
#if __DEBUG__
PrintMatrixOfSynapticWeights();
#endif
// = = = = = = = = = = = = = = = = = = = = = =
// = DELAY MATRIX
// = = = = = = = = = = = = = = = = = = = = = =
#if(GENERATE_NETWORK_FROM_EXTERNAL_DATA)
ImportDelayMatrixFromFile(INFILE_DELAY_MATRIX);
#else
for( int n = 0; n < NUM_TOTAL_NEURONS; ++n ) {
short synapseCorrespondingToDelay = 0;
for( int delayIdx = 0; delayIdx < MAX_SYNAPSE_DELAY; ++delayIdx ) {
numEntriesPerDelay[n][delayIdx] = 0;
}
if( n < NUM_EXCITATORY_NEURONS ) {
// For the excitatory neurons the delays are drawn from a uniform integer distribution.
// e.g.:
// delayIdx: 0 1 2 3 4 5 6 .......... each slot corresponds to a delay value, i.e., 0 -> 1 ms, 1 -> 2 ms etc.
// 0, 1 2, 3 4, 5 6, 7 8, 9 9, 10 11, 12 ... each slot has space for NUM_SYNAPSES_PER_NEURON (which is not needed) and
// holds the synapse numbers corresponding to that delay
// REFACTOR COMMENT: In the Matlab implementation this is random.
// This algorithm implicitly creates as many entries in delays as NUM_SYNAPSES_PER_NEURON.
// This only works if NUM_SYNAPSES_PER_NEURON is a multiple of MAX_SYNAPSE_DELAY!
for( int delayIdx = 0; delayIdx < MAX_SYNAPSE_DELAY; ++delayIdx ) {
for( int synListIdx = 0; synListIdx < (NUM_SYNAPSES_PER_NEURON / MAX_SYNAPSE_DELAY); ++synListIdx ) {
listOfSynsByNeuronAndDelay[n][delayIdx][synListIdx] = synapseCorrespondingToDelay++;
numEntriesPerDelay[n][delayIdx]++;
}
}
}
else {
// the inhibitory synapse delay is always 1 ms, thus, all synapses are added to delayIdx 0.
for( int synListIdx = 0; synListIdx < NUM_SYNAPSES_PER_NEURON; ++synListIdx ) {
listOfSynsByNeuronAndDelay[n][0][synListIdx] = synapseCorrespondingToDelay++;
numEntriesPerDelay[n][0]++;
}
}
}
#endif
// for all neurons: find the excitatory presynaptic neurons and delays and store them in lists
for( int n = 0; n < NUM_TOTAL_NEURONS; ++n ) {
numPreSynapticNeuronsOfTarget[n] = 0;
}
for( int n = 0; n < NUM_TOTAL_NEURONS; ++n ) {
for( int excPreSynNeuron = 0; excPreSynNeuron < NUM_EXCITATORY_NEURONS; ++excPreSynNeuron ) {
for( int s = 0; s < NUM_SYNAPSES_PER_NEURON; ++s ) {
int targetNeuron = matrixOfPostSynapticNeurons[excPreSynNeuron][s];
if( targetNeuron == n ) {
int idx = numPreSynapticNeuronsOfTarget[n];
listOfPresynapticNeurons[n][idx] = excPreSynNeuron;
// determine the delay of a connection and store it in a list
for( int delayIdx = 0; delayIdx < MAX_SYNAPSE_DELAY; ++delayIdx ) {
int excPreSynNeuron_numEntriesOfDelay = numEntriesPerDelay[excPreSynNeuron][delayIdx];
for( int i = 0; i < excPreSynNeuron_numEntriesOfDelay; ++i ) {
short synapse = listOfSynsByNeuronAndDelay[excPreSynNeuron][delayIdx][i];
int postSynNeuron = matrixOfPostSynapticNeurons[excPreSynNeuron][synapse];
if( postSynNeuron == n ) {
int idx = numPreSynapticNeuronsOfTarget[n];
listOfPresynapticDelays[n][idx] = delayIdx;
// REFACTOR COMMENT: delayIdx represents the delay value, where 0 represents 1ms, thus,
// one needs to distinguish between a 1ms delay and an empty entry.
// This is done by tracking the number of entries in delays_numEntriesPerDelay.
}
}
}
// maintain list of pointers into the weight matrices; used for weight update
{
int idx = numPreSynapticNeuronsOfTarget[n];
// REFACTOR COMMENT: dead code removed
// listOfPointersToSynapticWeights[n][idx] = &matrixOfSynapticWeights[excPreSynNeuron][s];
listOfPointersToSynapticWeights_derivatives[n][idx] = &matrixOfSynapticWeights_derivatives[excPreSynNeuron][s];
idx++;
numPreSynapticNeuronsOfTarget[n] = idx;
}
}
}
}
}
#if __DEBUG__
PrintMatrixOfSynapticDelays();
#endif
#ifdef OUTFILE_DELAYS
ExportDelayMatrixToFile( OUTFILE_DELAYS );
#endif
// REFACTOR COMMENT: the array is much larger and remains uninitialized
for( int n = 0; n < NUM_TOTAL_NEURONS; ++n ) {
for( int d = 0; d < 1 + MAX_SYNAPSE_DELAY; ++d ) {
LTP[n][d] = 0.0;
}
}
for( int n = 0; n < NUM_TOTAL_NEURONS; ++n ) {
LTD[n] = 0.0;
}
// REFACTOR COMMENT: just for the algorithm; does not contribute to the network activity
numFirings = 1; // dummy spike ...
firings[0][TIME] = -MAX_SYNAPSE_DELAY; // ... at -MAX_SYNAPSE_DELAY for ...
firings[0][NEURON] = 0; // ... neuron n = 0
#ifdef OUTFILE_CON_WEIGHT_DELAY_INITIAL
ExportConnectionMatrixWeightAndDelay(OUTFILE_CON_WEIGHT_DELAY_INITIAL);
#endif
}
// = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =
// = Initialize simulation: set seed, delete previously created files, open files
// = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =
void InitializeSimulation() {
srand( 0 ); // set a seed (repeatability)
#if( LOG_MIN_MAX_V_U )
v_max = 0.0;
v_min = 0.0;
u_max = 0.0;
u_min = 0.0;
#endif
#ifdef OUTFILE_STIMULUS
DeleteFile( OUTFILE_STIMULUS );
#endif
#ifdef OUTFILE_CONNECTIONS
DeleteFile( OUTFILE_CONNECTIONS );
#endif
#ifdef OUTFILE_DELAYS
DeleteFile( OUTFILE_DELAYS );
#endif
#ifdef OUTFILE_CON_WEIGHT_DELAY_INITIAL
DeleteFile(OUTFILE_CON_WEIGHT_DELAY_INITIAL);
#endif
#ifdef OUTFILE_WEIGHTS_INITIAL
DeleteFile(OUTFILE_WEIGHTS_INITIAL);
#endif
DeleteFile( OUTFILE_FIRINGS );
// open files which are required throughout the entire simulation
#ifdef OUTFILE_STIMULUS
pFileStimulusOutput = fopen( OUTFILE_STIMULUS, "w" );
if( pFileStimulusOutput == nullptr ) {
printf( "[ERROR] Failed to open stimulus output file. Filename: %s\n", OUTFILE_STIMULUS );
exit( RC_ERROR_EXIT );
}
#endif
#if(USE_EXTERNAL_STIMULUS)
pFileStimulusInput = fopen(INFILE_STIMULUS,"r" );
if(pFileStimulusInput == nullptr) {
printf( "[ERROR] Failed to open stimulus input file. Filename: %s\n", INFILE_STIMULUS );
exit( RC_ERROR_EXIT );
}
#endif
}
// = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =
// = Finalize simulation: free resources
// = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =
void FinalizeSimulation() {
if( pFileStimulusInput != nullptr ) fclose( pFileStimulusInput );
if( pFileStimulusOutput != nullptr ) fclose( pFileStimulusOutput );
}
#else // original version
// original version of the Izhikevich polychronization model available for download at
// https://www.izhikevich.org/publications/spnet.htm
#include "poly_spnet.cpp"
#endif