aeif_cond_alpha_multisynapse.h

Go to the documentation of this file.

/*
 *  aeif_cond_alpha_multisynapse.h
 *
 *  This file is part of NEST.
 *
 *  Copyright (C) 2004 The NEST Initiative
 *
 *  NEST 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.
 *
 *  NEST 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 NEST.  If not, see <http://www.gnu.org/licenses/>.
 *
 */
 
#ifndef AEIF_COND_ALPHA_MULTISYNAPSE_H
#define AEIF_COND_ALPHA_MULTISYNAPSE_H
 
// Generated includes:
#include "config.h"
#include <sstream>
 
#ifdef HAVE_GSL
 
// External includes:
#include <gsl/gsl_errno.h>
#include <gsl/gsl_matrix.h>
#include <gsl/gsl_odeiv.h>
 
// Includes from nestkernel:
#include "archiving_node.h"
#include "connection.h"
#include "event.h"
#include "nest_types.h"
#include "ring_buffer.h"
#include "universal_data_logger.h"
 
/* BeginUserDocs: neuron, integrate-and-fire, adaptation, conductance-based, soft threshold
 
Short description
+++++++++++++++++
 
Conductance based adaptive exponential integrate-and-fire neuron model
 
Description
+++++++++++
 
``aeif_cond_alpha_multisynapse`` is a conductance-based adaptive
exponential integrate-and-fire neuron model according to Brette and
Gerstner (2005) with multiple synaptic rise time and decay time
constants, and synaptic conductance modeled by an alpha function.
 
It allows an arbitrary number of synaptic time constants. Synaptic
conductance is modeled by an alpha function, as described by A. Roth
and M. C. W. van Rossum in Computational Modeling Methods for
Neuroscientists, MIT Press 2013, Chapter 6.
 
The time constants are supplied by an array, ``tau_syn``, and the pertaining
synaptic reversal potentials are supplied by the array ``E_rev``. Port numbers
are automatically assigned in the range from 1 to n_receptors.
During connection, the ports are selected with the property ``receptor_type``.
 
When connecting to conductance-based multisynapse models, all synaptic weights
must be non-negative.
 
The membrane potential is given by the following differential equation:
 
.. math::
 
 C dV/dt = -g_L(V-E_L) + g_L \cdot \Delta_T \cdot \exp((V-V_T)/\Delta_T)
 + I_{syn_{tot}}(V, t)- w + I_e
 
where
 
.. math::
 
 I_{syn_{tot}}(V,t) = \sum_i g_i(t) (V - E_{rev,i}) ,
 
the synapse `i` is excitatory or inhibitory depending on the value of
:math:`E_{rev,i}` and the differential equation for the
spike-adaptation current `w` is
 
.. math::
 
 \tau_w \cdot dw/dt = a(V - E_L) - w
 
When the neuron fires a spike, the adaptation current :math:`w <- w + b`.
 
For implementation details see the
`aeif_models_implementation <../model_details/aeif_models_implementation.ipynb>`_ notebook.
 
.. note::
 
    The default refractory period for ``aeif`` models is zero, consistent with the model definition in
    Brette & Gerstner :footcite:p:`Brette2005`.  Thus, an ``aeif`` neuron with default parameters can fire multiple
    spikes in a single time step, which can lead to exploding spike numbers and extreme slow-down of simulations.
 
    To avoid such unphysiological behavior, you should set a refractory time ``t_ref > 0``.
 
Parameters
++++++++++
 
The following parameters can be set in the status Dictionary.
 
======== ======= =======================================
**Dynamic state variables:**
--------------------------------------------------------
 V_m     mV      Membrane potential
 w       pA      Spike-adaptation current
======== ======= =======================================
 
======== ======= =======================================
**Membrane Parameters**
--------------------------------------------------------
 C_m     pF      Capacity of the membrane
 t_ref   ms      Duration of refractory period
 V_reset mV      Reset value for V_m after a spike
 E_L     mV      Leak reversal potential
 g_L     nS      Leak conductance
 I_e     pA      Constant external input current
 Delta_T mV      Slope factor
 V_th    mV      Spike initiation threshold
 V_peak  mV      Spike detection threshold
======== ======= =======================================
 
======== ======= ==================================
**Spike adaptation parameters**
---------------------------------------------------
 a       ns      Subthreshold adaptation
 b       pA      Spike-triggered adaptation
 tau_w   ms      Adaptation time constant
======== ======= ==================================
 
======== ============= ========================================================
**Synaptic parameters**
-------------------------------------------------------------------------------
E_rev    list of mV    Reversal potential
tau_syn  list of ms    Time constant of synaptic conductance
======== ============= ========================================================
 
============= ======= =========================================================
**Integration parameters**
-------------------------------------------------------------------------------
gsl_error_tol real    This parameter controls the admissible error of the
                      GSL integrator. Reduce it if NEST complains about
                      numerical instabilities.
============= ======= =========================================================
 
Sends
+++++
 
SpikeEvent
 
Receives
++++++++
 
SpikeEvent, CurrentEvent, DataLoggingRequest
 
References
++++++++++
 
.. footbibliography::
 
See also
++++++++
 
aeif_cond_alpha_multisynapse
 
Examples using this model
+++++++++++++++++++++++++
 
.. listexamples:: aeif_cond_alpha_multisynapse
 
EndUserDocs */
 
namespace nest
{
extern "C" int aeif_cond_alpha_multisynapse_dynamics( double, const double*, double*, void* );
 
void register_aeif_cond_alpha_multisynapse( const std::string& name );
 
class aeif_cond_alpha_multisynapse : public ArchivingNode
{
 
public:
  aeif_cond_alpha_multisynapse();
  aeif_cond_alpha_multisynapse( const aeif_cond_alpha_multisynapse& );
  ~aeif_cond_alpha_multisynapse() override;
 
  friend int aeif_cond_alpha_multisynapse_dynamics( double, const double*, double*, void* );
 
  using Node::handle;
  using Node::handles_test_event;
 
  size_t send_test_event( Node&, size_t, synindex, bool ) override;
 
  void handle( SpikeEvent& ) override;
  void handle( CurrentEvent& ) override;
  void handle( DataLoggingRequest& ) override;
 
  size_t handles_test_event( SpikeEvent&, size_t ) override;
  size_t handles_test_event( CurrentEvent&, size_t ) override;
  size_t handles_test_event( DataLoggingRequest&, size_t ) override;
 
  void get_status( Dictionary& ) const override;
  void set_status( const Dictionary& ) override;
 
private:
  void init_buffers_() override;
  void pre_run_hook() override;
  void update( Time const&, const long, const long ) override;
 
  // The next three classes need to be friends to access the State_ class/member
  friend class DynamicRecordablesMap< aeif_cond_alpha_multisynapse >;
  friend class DynamicUniversalDataLogger< aeif_cond_alpha_multisynapse >;
  friend class DataAccessFunctor< aeif_cond_alpha_multisynapse >;
 
  // ----------------------------------------------------------------
 
  struct Parameters_
  {
    double V_peak_;   
    double V_reset_;  
    double t_ref_;    
 
    double g_L;      
    double C_m;      
    double E_L;      
    double Delta_T;  
    double tau_w;    
    double a;        
    double b;        
    double V_th;     
 
    std::vector< double > tau_syn;  
    std::vector< double > E_rev;    
 
    double I_e;  
 
    double gsl_error_tol;  
 
    // boolean flag which indicates whether the neuron has connections
    bool has_connections_;
 
    Parameters_();  
 
    void get( Dictionary& ) const;              
    void set( const Dictionary&, Node* node );  
 
    inline size_t
    n_receptors() const
    {
      return E_rev.size();
    }
  };
 
  // ----------------------------------------------------------------
 
  struct State_
  {
    enum StateVecElems
    {
      V_M = 0,
      W,   // 1
      DG,  // 2
      G,   // 3
      STATE_VECTOR_MIN_SIZE
    };
 
    static const size_t NUMBER_OF_FIXED_STATES_ELEMENTS = 2;  // V_M, W
    static const size_t NUM_STATE_ELEMENTS_PER_RECEPTOR = 2;  // DG, G
 
    std::vector< double > y_;  
    int r_;                    
 
    State_( const Parameters_& );  
 
    void get( Dictionary& ) const;
    void set( const Dictionary&, Node* node );
 
  };  // State_
 
  // ----------------------------------------------------------------
 
  struct Buffers_
  {
    Buffers_( aeif_cond_alpha_multisynapse& );
    Buffers_( const Buffers_&, aeif_cond_alpha_multisynapse& );
 
    DynamicUniversalDataLogger< aeif_cond_alpha_multisynapse > logger_;
 
    std::vector< RingBuffer > spikes_;
    RingBuffer currents_;
 
    gsl_odeiv_step* s_;     
    gsl_odeiv_control* c_;  
    gsl_odeiv_evolve* e_;   
    gsl_odeiv_system sys_;  
 
    // Since IntegrationStep_ is initialized with step_, and the resolution
    // cannot change after nodes have been created, it is safe to place both
    // here.
    double step_;             
    double IntegrationStep_;  
 
    double I_stim_;
  };
 
  // ----------------------------------------------------------------
 
  struct Variables_
  {
    std::vector< double > g0_;
 
    double V_peak;
 
    unsigned int refractory_counts_;
  };
 
  // Data members -----------------------------------------------------------
 
  Parameters_ P_;
  State_ S_;
  Variables_ V_;
  Buffers_ B_;
  // Access functions for UniversalDataLogger -------------------------------
 
  DynamicRecordablesMap< aeif_cond_alpha_multisynapse > recordablesMap_;
 
  // Data Access Functor getter
  DataAccessFunctor< aeif_cond_alpha_multisynapse > get_data_access_functor( size_t elem );
  inline double
  get_state_element( size_t elem )
  {
    return S_.y_[ elem ];
  };
 
  // Utility function that inserts the synaptic conductances to the
  // recordables map
 
  std::string get_g_receptor_name( size_t receptor );
  void insert_conductance_recordables( size_t first = 0 );
};
 
inline size_t
aeif_cond_alpha_multisynapse::send_test_event( Node& target, size_t receptor_type, synindex, bool )
{
  SpikeEvent e;
  e.set_sender( *this );
 
  return target.handles_test_event( e, receptor_type );
}
 
inline size_t
aeif_cond_alpha_multisynapse::handles_test_event( CurrentEvent&, size_t receptor_type )
{
  if ( receptor_type != 0 )
  {
    throw UnknownReceptorType( receptor_type, get_name() );
  }
  return 0;
}
 
inline size_t
aeif_cond_alpha_multisynapse::handles_test_event( DataLoggingRequest& dlr, size_t receptor_type )
{
  if ( receptor_type != 0 )
  {
    throw UnknownReceptorType( receptor_type, get_name() );
  }
  return B_.logger_.connect_logging_device( dlr, recordablesMap_ );
}
 
inline void
aeif_cond_alpha_multisynapse::get_status( Dictionary& d ) const
{
  P_.get( d );
  S_.get( d );
  ArchivingNode::get_status( d );
 
  d[ names::recordables ] = recordablesMap_.get_list();
}
 
}  // namespace
 
#endif  // HAVE_GSL
#endif  // AEIF_COND_ALPHA_MULTISYNAPSE_H //