;;    Destexhe, A., Mainen, Z.F. and Sejnowski, T.J.  An efficient
;;    method for computing synaptic conductances based on a kinetic
;;    model of receptor binding Neural Computation 6: 10-14, 1994.


;; 1. during the pulse (from t=t0 to t=t1), C = Cmax, which gives:
;;
;;   R(t-t0) = Rinf + [ R(t0) - Rinf ] * exp (- (t-t0) / Rtau )		
;;
;; where 
;;   Rinf = Alpha * Cmax / (Alpha * Cmax + Beta) 
;; and
;;   Rtau = 1 / (Alpha * Cmax + Beta)
;;
;; 2. after the pulse (t>t1), C = 0, and one can write:
;;
;;    R(t-t1) = R(t1) * exp (- Beta * (t-t1) )				


(use srfi-1 mathh signal-diagram signal-diagram-dynamics )


(define (Destexhe94:construct t R Alpha Beta C_max C_dur Isyn gsyn Esyn Vsyn )

  `((f . ,(let* (
		 ;; on-regime: during the pulse, t \in [t0,t1]
		 (on-regime         (ASSIGN (R (+ R_ss (* (- R0 R_ss) (exp (neg (/ (- t t0) R_tau))))))
					    (endrelease (> t t1))
					    (spike spike)
					    (t0 t0)
					    (t1 t1)
					    (R0 R0)
					    (R1 R1)))
					    

		 ;; transient between off and on regimes: set t0 and t1 to indicate pulse boundaries
		 (off-on  (ASSIGN (R R)
				  (endrelease endrelease)
				  (spike spike)
				  (t0 t)
				  (t1 (+ t C_dur))
				  (R0 R)
				  (R1 R1)))

		 ;; off-regime: after the pulse, t > t1
		 (off-regime        (ASSIGN (R (* R1 (exp (neg (* (- t t1) Beta)))))
					    (endrelease false)
					    (spike spike)
					    (t0 t0)
					    (t1 t1)
					    (R0 R0)
					    (R1 R1)))

				    
 
		 ;; transient between on and off regimes: 
		 (on-off  (ASSIGN (R1 R)
				  (R  R)
				  (endrelease endrelease)
				  (spike spike)
				  (t0 t0)
				  (t1 t1)
				  (R0 R0)))

		 (synaptic_current (ASSIGN (Isyn (* gsyn R (- Vsyn Esyn)))))


		 )

	    (SEQUENCE
	     (RTRANSITION (TRANSIENT off-regime off-on 'spike )
			  (TRANSIENT on-regime on-off  'endrelease )
			  'spike 'endrelease 'firingstate)
	     synaptic_current)
	     ))

    (initial .
	     (
	      (h       1e-2)
	      (t       ,t)
	      (t0      ,t)
	      (t1      ,t)
	      (R       ,R)
	      (R0      ,R)
	      (R1      ,R)
	      (Alpha   ,Alpha)
	      (Beta    ,Beta)
	      (C_dur   ,C_dur)
	      (R_ss    ,(/ (* C_max Alpha) (+ (* C_max Alpha) Beta)))
	      (R_tau   ,(/ 1 (+ (* Alpha C_max) Beta)))
	      (Isyn    ,Isyn)
	      (gsyn    ,gsyn)
	      (Esyn    ,Esyn)
	      (Vsyn    ,Vsyn)
	      (spike   #f)
	      (endrelease #f)
	      (firingstate #f)
	     ))

    ))

(define (dataflow model)
  (let* ((f         (alist-ref 'f model)) 
	 (initial   (alist-ref 'initial model))
	 (input     (map car initial)))
    (print "input = " input)
    (print "dataflow f: "  (dataflow (construct f) input))
    ))


(define (codegen name model #!key (language 'scheme) (pre #t) (post #t))
  (let* ((f         (alist-ref 'f model)) 
	 (initial   (alist-ref 'initial model)))
    (case language
      ((scheme) (codegen/scheme name (construct f) initial: initial pre: pre solver: 'rk3))
      ((ML)     (codegen/ML name (construct f) initial: initial pre: pre post: post solver: 'rk3)))
    ))

;; t R Alpha Beta C_max C_dur Isyn gsyn Esyn Vsyn 
(define Destexhe94 (Destexhe94:construct 0.0 0.0 0.4 0.2 10.0 10.0 0.0 0.02 -75.0 0.0))


(with-output-to-file "Destexhe94_solver.scm"
  (lambda () (codegen 'Destexhe94 Destexhe94  language: 'scheme)))

(with-output-to-file "Destexhe94_solver.sml"
  (lambda () (codegen 'Destexhe94 Destexhe94 language: 'ML)))