/* This file was generated by {{nemoVersionString}} on {{currentTimestamp}} */

#include "nest.h"
#include "event.h"
#include "archiving_node.h"
#include "ring_buffer.h"
#include "connection.h"
#include "universal_data_logger.h"
#include "recordables_map.h"

{% if (ODEmethod == "cvode") %}

#include <sundials/sundials_types.h> /* definition of type realtype */
#include <nvector/nvector_serial.h>  /* serial N_Vector types, fcts., macros */
#include <cvode/cvode.h>             /* prototypes for CVODE fcts., consts. */
#include <cvode/cvode_diag.h>        /* prototype for CVDiag */
#define Ith(v,i)    NV_Ith_S(v,i)    /* Ith component in a vector */

{% elif (ODEmethod == "ida") %}

#include <sundials/sundials_types.h> /* definition of type realtype */
#include <nvector/nvector_serial.h>  /* serial N_Vector types, fcts., macros */
#include <ida/ida.h>                 /* prototypes for IDA fcts., consts. */
#include <ida/ida_dense.h>
#define Ith(v,i)    NV_Ith_S(v,i)    /* Ith component in a vector */

{% elif (ODEmethod == "gsl") %}

#include <gsl/gsl_errno.h>
#include <gsl/gsl_matrix.h>
#include <gsl/gsl_sf_exp.h>
#include <gsl/gsl_odeiv2.h>
#define Ith(v,i)    (v[i])

{% endif %}

{% if (SSmethod == "kinsol") %}

#include <sundials/sundials_types.h> /* definition of type realtype */
#include <nvector/nvector_serial.h>  /* serial N_Vector types, fcts., macros */
#include <kinsol/kinsol.h>           /* prototypes for KINSOL fcts. */
#include <kinsol/kinsol_dense.h>     /* prototype for KINDense */

{% else %}

#include <gsl/gsl_vector.h>
#include <gsl/gsl_multiroots.h>

{% endif %}


namespace nest {

{% if (ODEmethod == "cvode") %}
{% include "NEST-cvode-prelude.tmpl" %}
{% elif (ODEmethod == "ida") %}
{% include "NEST-ida-prelude.tmpl" %}
{% endif %}

{% if (SSmethod == "kinsol") %}
{% include "NEST-kinsol-prelude.tmpl" %}
{% endif %}

{% if (ODEmethod == "cvode") %}
  extern "C" int {{modelName}}_dynamics (double, const N_Vector, N_Vector, void*);
  extern "C" int {{modelName}}_event (double, N_Vector, double *, void*);
{% elif (ODEmethod == "ida") %}
  extern "C" int {{modelName}}_residual (double, N_Vector, N_Vector, N_Vector, void*);
  extern "C" int {{modelName}}_event (double, N_Vector, N_Vector, double *, void*);
{% else %}
  extern "C" int {{modelName}}_dynamics (double, const double*, double*, void*);
{% endif %}


{% if (steadyStateSize > 0) %}
{% if (SSmethod == "kinsol") %}
  extern "C" int {{modelName}}_steadystate (N_Vector, N_Vector, void*);

{% else %}
  extern "C" int {{modelName}}_steadystate (const gsl_vector *, void *, gsl_vector *);

{% endif %}
{% endif %}

  class {{modelName}} : public Archiving_Node { 

    public:

    ~{{modelName}} ();
    {{modelName}} (const {{modelName}} &);
    {{modelName}} ();

    using Node::connect_sender;
    using Node::handle;

    port check_connection(Connection&, port);
    
    void handle(SpikeEvent &);
    void handle(CurrentEvent &);
    void handle(DataLoggingRequest &); 
    
    port connect_sender(SpikeEvent &, port);
    port connect_sender(CurrentEvent &, port);
    port connect_sender(DataLoggingRequest &, port);
    
    void get_status(DictionaryDatum &) const;
    void set_status(const DictionaryDatum &);
    
    void init_node_(const Node& proto);
    void init_state_(const Node& proto);
    void init_buffers_();
    void calibrate();
    
    void update(Time const &, const long_t, const long_t);

    /**
     * Minimal spike receptor type.
     * @note Start with 1 so we can forbid port 0 to avoid accidental
     *       creation of connections with no receptor type set.
     */
    static const port MIN_SPIKE_RECEPTOR = 1;

    /** 
     * Spike receptors.
     */
    enum SpikeSynapseTypes { 

{% for synapticEvent in synapticEventDefs %}
{% if (loop.first) %}
         {{synapticEvent.pscName}}_SPIKE_RECEPTOR=MIN_SPIKE_RECEPTOR,

{% else %}
         {{synapticEvent.pscName}}_SPIKE_RECEPTOR,

{% endif %}
{% endfor %}
      SUP_SPIKE_RECEPTOR

    };
  
{% if (ODEmethod == "cvode") %}
    friend int {{modelName}}_dynamics (double, const N_Vector, N_Vector, void*);
    friend int {{modelName}}_event (double, N_Vector, double *, void*);
{% elif (ODEmethod == "ida") %}
    friend int {{modelName}}_residual (double, N_Vector, N_Vector, N_Vector, void*);
    friend int {{modelName}}_event (double, N_Vector, N_Vector, double *, void*);
{% elif (ODEmethod == "gsl") %}
    friend int {{modelName}}_dynamics (double, const double*, double*, void*);
{% endif %}

  
{% if (steadyStateSize > 0) %}
{% if (SSmethod == "kinsol") %}
   friend  int {{modelName}}_steadystate (N_Vector, N_Vector, void*);
{% else %}
   friend  int {{modelName}}_steadystate (const gsl_vector *, void *, gsl_vector *);
{% endif %}
{% endif %}

{% for synapticEvent in synapticEventDefs %}
   int {{synapticEvent.pscName}}_transients (long_t lag);
{% endfor %}  

   // The next two classes need to be friends to access the State_ class/member
   friend class RecordablesMap<{{modelName}}>;
   friend class UniversalDataLogger<{{modelName}}>;

  
   struct Parameters_ { 

{% if (not (parameters.parameterDefs == [])) %}
    double {% for p in parameters.parameterDefs %}{{p.name}}{% if (loop.last) %};{% else %},{% endif %}{% endfor %}
{% endif %}
{% if (not (parameters.defaultDefs == [])) %}
    double {% for p in parameters.defaultDefs %}{{p.name}}{% if (loop.last) %};{% else %},{% endif %}{% endfor %}
{% endif %}

    Parameters_();
    void get(DictionaryDatum&) const;
    void set(const DictionaryDatum&);

  }; // end struct Parameters_


  struct State_ { 

      enum StateVecElems {
{% for k in stateDefs %}
      {{k.name | upper}} = {{stateIndexMap[k.name]}}{% if (not (loop.last)) %},{% endif %}

{% endfor %}
      };

      double y_[{{stateSize}}]; 

      State_(const Parameters_& p); 
      State_(const State_& s);
      State_& operator=(const State_& s);
      void get(DictionaryDatum&) const;
      void set(const DictionaryDatum&, const Parameters_&);

{% if (ODEmethod == "gsl") %}
     int_t r_; /* refractory counter */
{% endif %}

  }; // end struct State_

{% if ((ODEmethod == "cvode") or (ODEmethod == "ida"))  %}
  struct Variables_ {};
{% elif (ODEmethod == "gsl")  %}
  struct Variables_ { int_t RefractoryCounts_; double U_old_; /* for spike-detection */ };
{% endif %}

  struct Buffers_ {

    Buffers_({{modelName}}&);
    Buffers_(const Buffers_&, {{modelName}}&);
    UniversalDataLogger<{{modelName}}> logger_;

{% if (ODEmethod == "cvode") %}
    N_Vector y;    //!< current state vector used by CVode
    void *   sys_;  //!< CVode control structure

{% elif (ODEmethod == "ida") %}
    N_Vector y, y1;    //!< current state vector used by IDA
    N_Vector yp;    //!< derivatives vector used by IDA
    void *   sys_;  //!< IDA control structure

{% elif (ODEmethod == "gsl") %}
    gsl_odeiv2_step*    s_;    //!< stepping function
    gsl_odeiv2_control* c_;    //!< adaptive stepsize control function
    gsl_odeiv2_evolve*  e_;    //!< evolution function
    gsl_odeiv2_system   sys_;  //!< struct describing system
    unsigned int N;  // size of state vector used by Jacobian
    double *u, *jac;  // intermediate state vectors used for Jacobian approximation

{% endif %}

{% for synapticEvent in synapticEventDefs %}
    RingBuffer spike_{{synapticEvent.pscId}};
{% endfor %}

    RingBuffer currents_;

    double_t step_;           //!< step size in ms
    double   IntegrationStep_;//!< current integration time step, updated by solver

  /** 
   * Input current injected by CurrentEvent.
   * This variable is used to transport the current applied into the
   * _dynamics function computing the derivative of the state vector.
   * It must be a part of Buffers_, since it is initialized once before
   * the first simulation, but not modified before later Simulate calls.
   */
    double_t I_stim_;

    }; // end struct Buffers_

  template <State_::StateVecElems elem>
  double_t get_y_elem_() const { return S_.y_[elem]; }

  Parameters_ P_;
  State_      S_;
  Variables_  V_;
  Buffers_    B_;

  static RecordablesMap<{{modelName}}> recordablesMap_;

  }; // end class {{modelName}}


  inline port {{modelName}}::check_connection(Connection& c, port receptor_type)
  {
    SpikeEvent e;
    e.set_sender(*this);
    c.check_event(e);
    return c.get_target()->connect_sender(e, receptor_type);
  }


  inline port {{modelName}}::connect_sender(SpikeEvent&, port receptor_type)
  {
    if ( receptor_type < MIN_SPIKE_RECEPTOR || receptor_type >= SUP_SPIKE_RECEPTOR )
    {
      if ( receptor_type < 0 || receptor_type >= SUP_SPIKE_RECEPTOR )
	throw UnknownReceptorType(receptor_type, get_name());
      else
	throw IncompatibleReceptorType(receptor_type, get_name(), "SpikeEvent");
    }
    return receptor_type;
  }
 
  inline port {{modelName}}::connect_sender(CurrentEvent&, port receptor_type)
  {
    if (receptor_type != 0)
      throw UnknownReceptorType(receptor_type, get_name());
    return 0;
  }

  inline port {{modelName}}::connect_sender(DataLoggingRequest& dlr, 
				            port receptor_type)
  {
    if (receptor_type != 0)
      throw UnknownReceptorType(receptor_type, get_name());
    return B_.logger_.connect_logging_device(dlr, recordablesMap_);
  }

  inline void {{modelName}}::get_status(DictionaryDatum &d) const
  {
    P_.get(d);
    S_.get(d);
    Archiving_Node::get_status(d);

    (*d)[names::recordables] = recordablesMap_.get_list();

    def<double_t>(d, names::t_spike, get_spiketime_ms());

    DictionaryDatum receptor_dict_ = new Dictionary();

{% for synapticEvent in synapticEventDefs %}
   (*receptor_dict_)[Name("{{synapticEvent.pscName}}")]  = {{synapticEvent.pscName}}_SPIKE_RECEPTOR;

{% endfor %}
               
    (*d)[names::receptor_types] = receptor_dict_;
  }


  inline void {{modelName}}::set_status(const DictionaryDatum &d)
  {
    Parameters_ ptmp = P_;  // temporary copy in case of errors
    ptmp.set(d);                       // throws if BadProperty
    State_      stmp = S_;  // temporary copy in case of errors
    stmp.set(d, ptmp);                 // throws if BadProperty

    // We now know that (ptmp, stmp) are consistent. We do not 
    // write them back to (P_, S_) before we are also sure that 
    // the properties to be set in the parent class are internally 
    // consistent.
    Archiving_Node::set_status(d);

    // if we get here, temporaries contain consistent set of properties
    P_ = ptmp;
    S_ = stmp;

    calibrate();
  }

} // end namespace nest