aeif_cond_beta_multisynapse.h

Go to the documentation of this file.

/*
 *  aeif_cond_beta_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_BETA_MULTISYNAPSE_H
#define AEIF_COND_BETA_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"
 
namespace nest
{
extern "C" int aeif_cond_beta_multisynapse_dynamics( double, const double*, double*, void* );
 
/* BeginUserDocs: neuron, adaptation, integrate-and-fire, conductance-based, soft threshold
 
Short description
+++++++++++++++++
 
Conductance based adaptive exponential integrate-and-fire neuron model
 
Description
+++++++++++
 
``aeif_cond_beta_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 a beta function.
 
It allows an arbitrary number of synaptic
rise time and decay time constants. Synaptic conductance is modeled by a
beta 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 two arrays, ``tau_rise`` and ``tau_decay`` for
the synaptic rise time and decay time, respectively. The synaptic
reversal potentials are supplied by the array ``E_rev``. The 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 `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_beta_multisynapse
 
EndUserDocs */
 
void register_aeif_cond_beta_multisynapse( const std::string& name );
 
class aeif_cond_beta_multisynapse : public ArchivingNode
{
 
public:
  aeif_cond_beta_multisynapse();
  aeif_cond_beta_multisynapse( const aeif_cond_beta_multisynapse& );
  ~aeif_cond_beta_multisynapse() override;
 
  friend int aeif_cond_beta_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_beta_multisynapse >;
  friend class DynamicUniversalDataLogger< aeif_cond_beta_multisynapse >;
  friend class DataAccessFunctor< aeif_cond_beta_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_rise;   
    std::vector< double > tau_decay;  
    std::vector< double > E_rev;      
 
    double I_e;  
 
    double gsl_error_tol;  
 
    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_beta_multisynapse& );
    Buffers_( const Buffers_&, aeif_cond_beta_multisynapse& );
 
    DynamicUniversalDataLogger< aeif_cond_beta_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_beta_multisynapse > recordablesMap_;
 
  // Data Access Functor getter
  DataAccessFunctor< aeif_cond_beta_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_beta_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_beta_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_beta_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_beta_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_BETA_MULTISYNAPSE_H //