;; 
;;
;; nemo documentation for the Chicken Scheme module system.
;;
;;
;; This program 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 3 of the
;; License, or (at your option) any later version.
;;
;; This program 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.
;;
;; A full copy of the GPL license can be found at
;; <http://www.gnu.org/licenses/>.
;;
;;


(require-library eggdoc srfi-1)
(import eggdoc srfi-1)

(define doc
  `((eggdoc:begin

     (name "nemo")

     (description "Neuron model description language.")

     (author (url "http://chicken.wiki.br/users/ivan-raikov" "Ivan Raikov"))

     (history
      (version "2.4"   "Converted to using getopt-long")
      (version "2.3"   "Added eggdoc as a dependency")
      (version "2.2"   "Added stx-engine.scm to file manifest")
      (version "2.1"   "Ported to Chicken 4")
      (version "2.0"   "Introduced functors")
      (version "1.15"  "Added nmodl-depend option")
      (version "1.14"  "Added support for exponential Euler integration")
      (version "1.13"  "Change in the integration method used for the AKP example")
      (version "1.12"  "Added support for binary conductances and conservation equations")
      (version "1.11"  "Bug fixes in the current equations part of NMODL code generator")
      (version "1.10"  "AKP06 example is now installed in CHICKEN-HOME/nemo/examples")
      (version "1.9"   "Documentation and example updates")
      (version "1.8"   "Bug fixes related to kinetic equation processing")
      (version "1.6"   "Added infix expression parser (nemo format)")
      (version "1.0"   "Initial release"))


     (requires 
	       (url "datatype.html"        "datatype")
	       (url "digraph.html"         "digraph")
	       (url "eggdoc.html"          "eggdoc")
	       (url "environments.html"    "environments")
	       (url "getopt-long.html"     "getopt-long")
	       (url "graph-bfs.html"       "graph-bfs")
	       (url "graph-cycles.html"    "graph-cycles")
	       (url "lalr.html"            "lalr")
	       (url "matchable.html"       "matchable")
	       (url "mathh.html"           "mathh")
	       (url "strictly-pretty.html" "strictly-pretty")
	       (url "ssax.html"            "ssax")
	       (url "sxpath.html"          "sxpath")
	       (url "sxml-transforms.html" "sxml-transforms")
	       (url "varsubst.html"        "varsubst")
	       (url "vector-lib.html"      "vector-lib")
	       )


     (usage "nemo [options...] [input files ...]")
     (download "nemo.egg")

     (documentation
      
      (p (tt "NEMO") " is a program that reads an on channel description "
	 "in a format based on ChannelML and generates a corresponding description "
	 "in the " (url "http://www.neuron.yale.edu/neuron/docs/help/neuron/nmodl/nmodl.html" "NMODL")
	 "language used by the " (url "http://www.neuron.yale.edu/neuron/" "NEURON simulator") ". ")

      (subsection "Options"
        (p (symbol-table
	    (describe "-i FORMAT"          "specify input format (nemo, xml, sxml, s-exp)")
	    (describe "--xml[=FILE]"       "write XML output to file (default: <model-name>.xml")
	    (describe "--sxml[=FILE]"      "write SXML output to file (default: <model-name>.sxml")
	    (describe "--nmodl[=FILE]"     "write NMODL output to file (default: <model-name>.mod")
	    (describe "--nmodl-method=METHOD"  "specify NMODL integration method (cnexp, derivimplicit)")
	    (describe "--nmodl-kinetic=[STATES]" "use NMODL kinetic equations for the given reactions")
	    (describe "--nmodl-depend=VARS" "specify DEPEND variables for NMODL interpolation tables ")
	    (describe "-t"                 "use interpolation tables in generated code")
	    (describe "-h, --help"         "print help"))))

      (subsection "Model description language"
	(p "The following constructs comprise the model description language: ")
        (p (ebnf-def "MODEL"
	    (ebnf-choice
	     ( (ebnf-form (ebnf-kw "INPUT ") 
			  (ebnf-choice (ebnf-var "ID") 
			     (ebnf-form (ebnf-var "ID") (ebnf-opt (ebnf-kw "AS ") (ebnf-var "LOCAL-ID")) 
					(ebnf-opt (ebnf-kw "FROM ") (ebnf-var "NAMESPACE"))))
			  "... )")
	       (comment "Declares one or several imported quantities. "
			"If the optional " (ebnf-kw "AS") " parameter is given, "
			"then the quantity is imported as " (ebnf-var "LOCAL-ID") ". "
			"If the optional " (ebnf-kw "FROM") " parameter is given, "
			"then the quantity is imported from namespace " (ebnf-var "NAMESPACE") ". "))
	       
	     ( (ebnf-form (ebnf-kw "OUTPUT") (ebnf-var "ID") )
	       (comment "Declares that an existing quantity be exported."))

	     ( (ebnf-form (ebnf-kw "CONST") (ebnf-var "ID") "=" (ebnf-var "EXPR"))
	       (comment "Declares a constant quantity (its value will be computed at declaration time)."))

	     ( (ebnf-form (ebnf-kw "DEFUN") (ebnf-var "ID") (ebnf-form (ebnf-var "ARG-ID") "..." )
			  (ebnf-var "EXPR"))
	       (comment "Declares a function (a parameterized expression with no free variables)."))

	     ( (ebnf-form (ebnf-var "ID") "=" (ebnf-var "EXPR"))
	       (comment "Declares an assigned quantity "
			"(an expression that can refer to other quantities in the system)."))

	     ( (ebnf-form (ebnf-kw "REACTION") (ebnf-var "ID") (ebnf-var "TRANSITIONS")
			  (ebnf-var "INITIAL-EXPR") (ebnf-var "OPEN-ID"))
	       (p "Declares a reaction quantity. See below for the syntax of state transition equations. "
		  (ebnf-var "INITIAL-EXPR") " is an expression that computes the initial value. "
		  (ebnf-var "OPEN-ID") " is the name of the open state. It must be one of the states "
		  "defined by the transition equations. "))

	     ( (ebnf-form (ebnf-kw "COMPONENT") 
			  (ebnf-form (ebnf-kw "TYPE") (ebnf-var "ID"))
			  (ebnf-form (ebnf-kw "NAME") (ebnf-var "ID"))
			  (ebnf-var "ELEMENTS"))
	       (p "Declares a system component (a quantity that can contain other quantities)."))
	     )))

	(subsubsection "Expressions"
          (p "Expressions in the model description language are defined as: ")
	  (p (ebnf-def "EXPR"
	       (ebnf-choice
		       ( (ebnf-var "NUM")
			 (comment "A numeric constant."))
		       ( (ebnf-var "ID")
			 (comment "A variable name."))
		       ( (ebnf-form (ebnf-var "ID") (ebnf-form (ebnf-var "EXPR") "..."))
			 (comment "A function invocation."))
		       ( (ebnf-form (ebnf-var "EXPR") (ebnf-var "OP") (ebnf-var "EXPR") )
			 (comment "Arithmetic operator invocation. The following operators are supported: "
				  (tt "+ - / * > < <= >= ^")))
		       ( (ebnf-form (ebnf-kw "LET") (ebnf-form (ebnf-var "BINDINGS")) (ebnf-var "EXPR") )
			 (comment "Local variables declaration. Each element in " 
				  (ebnf-var "BINDINGS") " is of the form: "
				  (ebnf-form (ebnf-var "ID") (ebnf-var "EXPR"))))
		       ( (ebnf-form (ebnf-kw "IF") (ebnf-var "CONDITION")
				    (ebnf-kw "THEN") (ebnf-var "EXPR") 
				    (ebnf-kw "ELSE") (ebnf-var "EXPR") )
			 (comment "Conditional expression. The expression after " (ebnf-kw "IF") 
				  " must be a comparison expression. ")))
		       )))

	(subsubsection "State transition equations"
          (p "State transition equations in the model description language are defined as: ")
	  (p (ebnf-def "TRANSITION"
               (ebnf-choice
		( (ebnf-form (ebnf-kw "->") (ebnf-var "SRC-ID") (ebnf-var "DEST-ID")
			     (ebnf-var "EXPR") )
		  (comment "Declares that a transition occurs from state "
			   (ebnf-var "SRC-ID") " to state " (ebnf-var "DEST-ID") 
			   " at rate computed by " (ebnf-var "EXPR") ". "))
		( (ebnf-form (ebnf-kw "<->") (ebnf-var "SRC-ID") (ebnf-var "DEST-ID")
			     (ebnf-var "EXPR-1") (ebnf-var "EXPR-2") )
		  (comment "Declares that a transition occurs from state "
			   (ebnf-var "SRC-ID") " to state " (ebnf-var "DEST-ID") 
			   " and vice versa, at rates computed by " (ebnf-var "EXPR-1")
			   " and " (ebnf-var "EXPR-2") ". "))
		))))

	(subsubsection "Ion channel definitions"
          (p "Currently, the " (tt "NMODL") " code generator recognizes and generates code for "
	     "ion channel components that are defined as follows: ")
	  (p (ebnf-form (ebnf-compound "COMPONENT (TYPE ion-channel) (NAME {NAME})"

  	       ( (ebnf-form (ebnf-kw "COMPONENT") (ebnf-form (ebnf-kw "TYPE") "gate") "..." )
		 (comment "One or more gate definitions. Each component of type gate "
			  "must export the reactions that characterize the gate dynamics. "))
	       
	       ( (ebnf-form (ebnf-kw "COMPONENT") 
			    (ebnf-form (ebnf-kw "TYPE") "pore") "...")
		 (comment "Conductance law definition. This component must export a "
			  "constant maximal conductance, or an assigned quantity whose "
			  "equation represents the conductance law used. "))
	       
	       ( (ebnf-opt (ebnf-form (ebnf-kw "COMPONENT") (ebnf-form (ebnf-kw "TYPE") "permeating-substance") "..."))
		 (comment ""))

	       ( (ebnf-opt (ebnf-form (ebnf-kw "COMPONENT") (ebnf-form (ebnf-kw "TYPE") "accumulating-substance") "..."))
		 (comment ""))
			    
	       ))))

	(subsubsection "Hodgkin-Huxley ionic conductance extension"
	  (p "The Hodgkin-Huxley ionic conductance extension is a shortcut that declares "
	     "a reaction corresponding to the Hodgkin-Huxley formulation "
	     "of ion channel dynamics. ")

	  (p (ebnf-form (ebnf-compound "HH-IONIC-GATE"
  	       ( (ebnf-form (ebnf-var "ION-NAME" )
			    (comment "Ion name: exported variables will be of the form " (tt "{ion}_{id}") ". ")
			    
			    (ebnf-form (ebnf-kw "M-POWER") (ebnf-var "INTEGER") )
			    (comment "Power of state variable " (tt "M") ". " )

			    (ebnf-form (ebnf-kw "H-POWER") (ebnf-var "INTEGER") )
			    (comment "Power of state variable " (tt "H") ". If zero, the initial value and "
				     "equations for this variable can be omitted. ")

			    (ebnf-form (ebnf-kw "INITIAL-M") (ebnf-var "EXPR") )
			    (comment "Expression that computes initial value for state variable " (tt "M") ". " )

			    (ebnf-form (ebnf-kw "INITIAL-H") (ebnf-var "EXPR") )
			    (comment "Expression that computes initial value for state variable " (tt "H") ". " )


			    (ebnf-form (ebnf-kw "M-ALPHA") (ebnf-var "EXPR") )
			    (comment "Closed state to open state rate expression for state variable " (tt "M") ". " )

			    (ebnf-form (ebnf-kw "M-BETA") (ebnf-var "EXPR") )
			    (comment "Open state to closed state rate expression for state variable " (tt "M") ". " )

			    (ebnf-form (ebnf-kw "H-ALPHA") (ebnf-var "EXPR") )
			    (comment "Closed state to open state rate expression for state variable " (tt "H") ". " )

			    (ebnf-form (ebnf-kw "H-BETA") (ebnf-var "EXPR") )
			    (comment "Open state to closed state rate expression for state variable " (tt "H") ". " )


			    (ebnf-form (ebnf-kw "M-INF") (ebnf-var "EXPR") )
			    (comment "Steady state expression for variable " (tt "M") ". " )

			    (ebnf-form (ebnf-kw "M-TAU") (ebnf-var "EXPR") )
			    (comment "Time constant expression for variable " (tt "M") ". " )

			    (ebnf-form (ebnf-kw "H-INF") (ebnf-var "EXPR") )
			    (comment "Steady state expression for variable " (tt "H") ". " )

			    (ebnf-form (ebnf-kw "H-TAU") (ebnf-var "EXPR") )
			    (comment "Time constant expression for variable " (tt "H") ". " )


			    ))
		 )))
	  )

	)

      )
     (examples (pre #<<EOF

;; Akemann and Knoepfel, J.Neurosci. 26 (2006) 4602
(nemo-model AKP06

  ((input v 
	  (cai from ion-pools)
	  (ica from ion-currents))
	     
   (const Vrest   = -68)
   (const diam    = 20)
   (const celsius = 24)
   
   (const temp_adj = (pow (3 ((celsius - 22) / 10))))


   (defun ghk (v celsius ci co)
     (let ((F 96485.0) (R 8.3145))
       (let ((zeta ((2e-3 * F * v) / (R * (273.19 + celsius)))))
	 (if ((abs (1.0 - exp (neg (zeta)))) < 1e-6)
	     then (1e-6 * (2 * F) * (ci - (co * exp (neg (zeta)))) * (1.0 + (zeta / 2)))
	     else ((1e-6 * (2 * zeta * F) * (ci - (co * exp (neg (zeta))))) / (1.0 - exp (neg (zeta))))))))

   (component (type membrane-capacitance)
	   (const C_m = 1)
	   (output C_m))

   (component (type ion-channel) (name Kv1)
	      
	      (component (type gate)
			 
			 ;; rate functions
			 (defun Kv1_amf (v) 
			   (let ((cma 0.12889)
				 (cka -33.90877)
				 (cva 45))
			     (cma * (exp (neg ((v + cva) / cka))))))
			 
			 (defun Kv1_bmf (v)
			   (let ((cmb   0.12889)
				 (ckb   12.42101)
				 (cvb   45))
			     (cmb * (exp (neg ((v + cvb) / ckb))))))
			 
			 (hh-ionic-gate 
			  (Kv1  ;; ion name: exported variables will be of the form {ion}_{id}
			   (initial-m (Kv1_amf (Vrest) / (Kv1_amf (Vrest) + Kv1_bmf (Vrest))))
			   (m-power   4)
			   (h-power   0)
			   (m-alpha   (temp_adj * Kv1_amf (v) ))
			   (m-beta    (temp_adj * Kv1_bmf (v) ))))
			 
			 )
	      
	      (component (type pore)
			 (const  gbar_Kv1  = 0.011)
			 (output gbar_Kv1 ))
	      
	      (component (type permeating-substance) (name k)
			 (const e_Kv1 = -85)
			 (output e_Kv1 ))
	      
	      ) ;; end Kv1 current
   
   (component (type ion-channel) (name Kv3)
	      
	      (component (type gate)
			 
			 ;; rate functions
			 (defun Kv3_amf (v) 
			   (let ((ca   0.22)
				 (cva  16)
				 (cka  -26.5))
			     (ca * exp((neg (v + cva)) / cka) )))
			 
			 (defun Kv3_bmf (v) 
			   (let ((ca   0.22)
				 (cvb  16)
				 (ckb  26.5))
			     (ca * exp((neg (v + cvb)) / ckb) )))
			 
			 (hh-ionic-gate 
			  (Kv3  ;; ion name: exported variables will be of the form {ion}_{id}
			   (initial-m (Kv3_amf (Vrest) / (Kv3_amf (Vrest) + Kv3_bmf (Vrest))))
			   (m-power   4)
			   (h-power   0)
			   (m-alpha   (temp_adj * Kv3_amf (v) ))
			   (m-beta    (temp_adj * Kv3_bmf (v) ))))
			 
			 )
	      
	      (component (type pore)
			 (const  gbar_Kv3  = 0.005)
			 (output gbar_Kv3 ))
	      
	      (component (type permeating-substance) (name k)
			 (const e_Kv3 = -85)
			 (output e_Kv3 ))

	      (component (type binary-gate) 
			 (const switch_Kv3 = 0)
			 
			 (const e0 = 1.60217646e-19)
			 (const gunit = 16)
			 (const nc = (1e12 * gbar_Kv3 / gunit))
			 (const zn = 1.9196)

			 (defun gate_flip_Kv3 (v m) 
			   (let ((a    (Kv3_amf (v)))
				 (ab   (a + Kv3_bmf (v)))
				 (tau  (1 / ab))
				 (inf  (a / ab)))
			     ((inf - m) / tau)))

			 (i_gate_Kv3 = 
                           (if (switch_Kv3 > 0) then (nc * 1e6 * e0 * 4 * zn * gate_flip_Kv3(v Kv3_m)) else 0))

			 (output i_gate_Kv3 ))

	      ) ;; end Kv3 current



   (component (type ion-channel) (name Kv4)

	      (component (type gate)
			 
			 ;; rate functions
			 
			 
			 (defun Kv4_amf (v)  
			   (let ((can    0.15743)
				 (ckan   -32.19976)
				 (cvan   57))
			     (can * exp (neg ((v + cvan) / ckan)))))
			 

			 (defun Kv4_bmf (v)  
			   (let ((cbn   0.15743)
				 (ckbn  37.51346)
				 (cvbn  57))
			     (cbn * exp (neg ((v + cvbn) / ckbn)))))

			 
			 (defun Kv4_ahf (v)  
			   (let ((cah   0.01342)
				 (ckah  -7.86476)
				 (cvah  60))
			     (cah / (1.0 + (exp (neg ((v + cvah) / ckah)))))))

			 
			 (defun Kv4_bhf (v)  
			   (let ((cbh   0.04477)
				 (ckbh  11.3615)
				 (cvbh  54))
			     (cbh / (1.0 + (exp (neg ((v + cvbh) / ckbh)))))))

			 (hh-ionic-gate 
			  (Kv4  ;; ion name: exported variables will be of the form {ion}_{id}
			   (initial-m (Kv4_amf (Vrest) / (Kv4_amf (Vrest) + Kv4_bmf (Vrest))) )
			   (initial-h (Kv4_ahf (Vrest) / (Kv4_ahf (Vrest) + Kv4_bhf (Vrest))) ) 
			   (m-power   4)
			   (h-power   1)
			   (m-alpha   (temp_adj * Kv4_amf (v)))
			   (m-beta    (temp_adj * Kv4_bmf (v)))
			   (h-alpha   (temp_adj * Kv4_ahf (v)))
			   (h-beta    (temp_adj * Kv4_bhf (v)))
			   ))
			 
			 )

	      (component (type pore)
			 (const  gbar_Kv4  = 0.0039)
			 (output gbar_Kv4 ))
	      
	      (component (type permeating-substance) (name k)
			 (const e_Kv4 = -85)
			 (output e_Kv4 ))
	      
	      ) ;; end Kv4 current


   (component (type ion-channel) (name Ih)
	      
	      (component (type gate)
			 
			 ;; rate functions
			 
			 (defun Ih_inf (v)
			   (let ((cvn   90.1)
				 (ckn   -9.9))
			     (1.0 / (1.0 + exp (neg ((v + cvn) / ckn) )))))
			 
			 (defun Ih_tau (v)   
			   (let ((cct  190)
				 (cat  720)
				 (cvt  81.5)
				 (ckt  11.9))
			     (cct + (cat * exp (neg (pow (((v + cvt) / ckt) 2)))))))
			 
			 (hh-ionic-gate 
			  (Ih  ;; ion name: exported variables will be of the form {ion}_{id}
			   (initial-m (Ih_inf (Vrest)))
			   (m-power   1)
			   (h-power   0)
			   (m-inf     (Ih_inf (v)))
			   (m-tau     (Ih_tau (v) / temp_adj))
			   ))
			 
			 )

	      (component (type pore)
			 (const  gbar_Ih  = 0.0002)
			 (output gbar_Ih ))
	      
	      (component (type permeating-substance) (name non-specific)
			 (const e_Ih = -30)
			 (output e_Ih ))
	      
	      ) ;; end Ih current

   (component (type ion-channel) (name Leak)
	      
	      (component (type pore)
			 (const  gbar_Leak  = 9e-5)
			 (output gbar_Leak ))
	      
	      (component (type permeating-substance) (name non-specific)
			 (const e_Leak = -61)
			 (output e_Leak ))
	      
	      ) ;; end leak current


   (component (type ion-channel) (name CaP)
	      
	      (component (type gate)
			 
			 ;; rate functions
			 (defun CaP_inf (v)  
			   (let ((cv  19)  (ck  5.5))
			     (1.0 / (1.0 + exp (neg ((v + cv) / ck))))))
			 
			 (defun CaP_tau (v)  
			   (if (v > -50) 
			       then (1e3 * (0.000191 + (0.00376 * pow ((exp (neg ((v + 41.9) / 27.8)))  2))))
			       else (1e3 * (0.00026367 + (0.1278 * exp (0.10327 * v))))))

			 (hh-ionic-gate 
			  (CaP  ;; ion name: exported variables will be of the form {ion}_{id}
			   (initial-m  (CaP_inf (Vrest)))
			   (m-power    1)
			   (h-power    0)
			   (m-inf      (CaP_inf (v)))
			   (m-tau      (CaP_tau (v) / temp_adj))))
			 
			 )

	      (component (type pore)
			 (const gmax_CaP  = 0.01667)
			 (const cao       = 2.4)
			 (gbar_CaP        = (gmax_CaP * ghk (v celsius cai cao)))
			 (output gbar_CaP ))
	      
	      (component (type accumulating-substance) (name ca) )
	      
	      
	      ) ;; end CaP current

   (component (type ion-channel) (name CaBK)
	      
	      (component (type gate)
			 
			 ;; rate functions
			 
			 (defun CaBK_zinf (ca)  
			   (let ((zhalf 0.001))
			     (1 / (1 + (zhalf / ca)))))

			 (const CaBK_ztau = 1.0)
			 
			 (defun CaBK_minf (v)
			   (let ((cvm  28.9)
				 (ckm  6.2))
			     (1.0 / (1.0 + exp (neg ((v + 5.0 + cvm) / ckm))))))
			 
			 (defun CaBK_mtau (v)
			   (let ((ctm     0.000505)
				 (cvtm1   86.4)
				 (cktm1   -10.1)
				 (cvtm2   -33.3)
				 (cktm2   10))
			     (ctm + (1.0 / (exp (neg ((v + 5.0 + cvtm1) / cktm1)) +
						exp (neg ((v + 5.0 + cvtm2) / cktm2)))))))
			 
			 (defun CaBK_hinf (v) 
			   (let ((ch    0.085)
				 (cvh   32)
				 (ckh   -5.8))
			     (ch + ((1.0 - ch) / (1.0 + (exp (neg ((v + 5.0 + cvh) / ckh))))))))
			 
			 (defun CaBK_htau (v)
			   (let ((cth     0.0019)
				 (cvth1   48.5)
				 (ckth1   -5.2)
				 (cvth2   -54.2)
				 (ckth2   12.9))
			     (cth + (1.0 / (exp (neg ((v + cvth1) / ckth1)) +
						exp (neg ((v + cvth2) / ckth2)))))))
			 
			 (reaction
			  (CaBK_z
			   (transitions (<-> C O (CaBK_zinf (cai) / CaBK_ztau)
					     ((1 - CaBK_zinf (cai)) / CaBK_ztau)))
			   (initial   (CaBK_zinf (1e-4)))
			   (open O)  (power 2)))
			 
			 (output CaBK_z )  

			 
			 (hh-ionic-gate 
			  (CaBK  ;; ion name: exported variables will be of the form {ion}_{id}
			   (initial-m  (CaBK_minf (Vrest) / temp_adj))
			   (initial-h  (CaBK_hinf (Vrest) / temp_adj))
			   (m-power    3)
			   (h-power    1)
			   (m-inf      (CaBK_minf (v) / temp_adj) )
			   (m-tau      (CaBK_mtau (v) / temp_adj) )
			   (h-inf      (CaBK_hinf (v) / temp_adj) )
			   (h-tau      (CaBK_htau (v) / temp_adj) )))
			 
			 )
	      
	      (component (type pore)
			 (const  gbar_CaBK  = 0.014)
			 (output gbar_CaBK ))
	      
	      (component (type permeating-substance) (name k)
			 (const e_CaBK = -85)
			 (output e_CaBK ))
	      
	      ) ;; end BK current


    (component (type decaying-pool) (name ca)
 	      (const  F = 96485.0)
 	      (const  ca_initial = 1e-4)
 	      (const  ca_depth   = 0.1)
 	      (const  ca_beta    = 1.0)
	      
 	      (d (ca) = ((neg (ica) / (2 * F * ca_initial * ca_depth)) -
 			 (ca * ca_beta * temp_adj))
 			(initial ca_initial))
	      
 	      (output ca)
 	      )
	    

   (functor (name Nafun) (type ion-channel) 
	    ( Na_Con Na_Coff Na_Oon Na_Ooff Na_alfac Na_btfac
		     Na_alpha Na_beta Na_gamma Na_delta Na_epsilon Na_zeta Na_x1 Na_x2
		     Na_x3 Na_x4 Na_x5 Na_x6 ) =
		     
		     (component (type gate)
				
		;; rate functions
		(f01 = (4.0 * Na_alpha * exp (v / Na_x1) * temp_adj))
		(f02 = (3.0 * Na_alpha * exp (v / Na_x1) * temp_adj))
		(f03 = (2.0 * Na_alpha * exp (v / Na_x1) * temp_adj))
		(f04 = (Na_alpha * exp (v / Na_x1) * temp_adj))
		(f0O = (Na_gamma * exp (v / Na_x3) * temp_adj))
		(fip = (Na_epsilon * exp (v / Na_x5) * temp_adj))
		(f11 = (4.0 * Na_alpha * Na_alfac * exp (v / Na_x1) * temp_adj))
		(f12 = (3.0 * Na_alpha * Na_alfac * exp (v / Na_x1) * temp_adj))
		(f13 = (2.0 * Na_alpha * Na_alfac * exp (v / Na_x1) * temp_adj))
		(f14 = (Na_alpha * Na_alfac * exp (v / Na_x1) * temp_adj))
		(f1n = (Na_gamma * exp (v / Na_x3) * temp_adj))
		
		(fi1 = (Na_Con * temp_adj))
		(fi2 = (Na_Con * Na_alfac * temp_adj))
		(fi3 = (Na_Con * Na_alfac * Na_alfac * temp_adj))
		(fi4 = (Na_Con * Na_alfac * Na_alfac * Na_alfac * temp_adj))
		(fi5 = (Na_Con * Na_alfac * Na_alfac * Na_alfac * Na_alfac * temp_adj))
		(fin = (Na_Oon * temp_adj))
		
		(b01 = (Na_beta * exp (v / Na_x2) * temp_adj))
		(b02 = (2.0 * Na_beta * exp (v / Na_x2) * temp_adj))
		(b03 = (3.0 * Na_beta * exp (v / Na_x2) * temp_adj))
		(b04 = (4.0 * Na_beta * exp (v / Na_x2) * temp_adj))
		(b0O = (Na_delta * exp (v / Na_x4) * temp_adj))
		(bip = (Na_zeta * exp (v / Na_x6) * temp_adj))
		
		(b11 = (Na_beta * Na_btfac * exp (v / Na_x2) * temp_adj))
		(b12 = (2.0 * Na_beta * Na_btfac * exp (v / Na_x2) * temp_adj))
		(b13 = (3.0 * Na_beta * Na_btfac * exp (v / Na_x2) * temp_adj))
		(b14 = (4.0 * Na_beta * Na_btfac * exp (v / Na_x2) * temp_adj))
		(b1n = (Na_delta * exp (v / Na_x4) * temp_adj))
		
		(bi1 = (Na_Coff * temp_adj))
		(bi2 = (Na_Coff * Na_btfac * temp_adj))
		(bi3 = (Na_Coff * Na_btfac * Na_btfac * temp_adj))
		(bi4 = (Na_Coff * Na_btfac * Na_btfac * Na_btfac * temp_adj))
		(bi5 = (Na_Coff * Na_btfac * Na_btfac * Na_btfac * Na_btfac * temp_adj))
		(bin = (Na_Ooff * temp_adj))
		
		(reaction
		    (Na_z
		     (transitions
		      (<-> C1 C2 f01 b01)
		      (<-> C2 C3 f02 b02)
		      (<-> C3 C4 f03 b03)
		      (<-> C4 C5 f04 b04)
		      (<-> C5 O  f0O b0O)
		      (<-> O  B  fip bip)
		      (<-> O  I6 fin bin)
		      (<-> I1 I2 f11 b11)
		      (<-> I2 I3 f12 b12)
		      (<-> I3 I4 f13 b13)
		      (<-> I4 I5 f14 b14)
		      (<-> I5 I6 f1n b1n)
		      (<-> C1 I1 fi1 bi1)
		      (<-> C2 I2 fi2 bi2)
		      (<-> C3 I3 fi3 bi3)
		      (<-> C4 I4 fi4 bi4)
		      (<-> C5 I5 fi5 bi5))

		      (conserve (1 = (C1 + C2 + C3 + C4 + C5 + O + B + I1 + I2 + I3 + I4 + I5 + I6)))
		     
		     (open O)   (power 1)))
		   
		   (output Na_z )  
		   
		   )
			
		(component (type pore)
			   (const  gbar  = 0.016)
			   (output gbar ))
		
		(component (type permeating-substance) (name na)
			   (const e = -88)
			   (output e ))
		
		) ;; end Nafun functor



     (component (name Na) =
		Nafun ((const Na_Con   = 0.005)
		       (const Na_Coff  = 0.5)
		       (const Na_Oon   = 2.3)
		       (const Na_Ooff  = 0.005)
		       
		       (const Na_alfac = (pow ((Na_Oon / Na_Con) (1.0 / 4.0))))
		       (const Na_btfac = (pow ((Na_Ooff / Na_Coff) (1.0 / 4.0))))
		       
		       (const Na_alpha = 150)
		       (const Na_beta  = 3)
		       (const Na_gamma = 150)
		       (const Na_delta = 40)
		       (const Na_epsilon = 1e-12)
		       (const Na_zeta = 0.03)
		       (const Na_x1 = 20)
		       (const Na_x2 = -20)
		       (const Na_x3 = 1000000000000.0)
		       (const Na_x4 = -1000000000000.0)
		       (const Na_x5 = 1000000000000.0)
		       (const Na_x6 = -25)))
     ;; end Na current

     (component (name Narsg) =
		Nafun ((const Na_Con   = 0.005)
		       (const Na_Coff  = 0.5)
		       (const Na_Oon   = 0.75)
		       (const Na_Ooff  = 0.005)
		       
		       (const Na_alfac = (pow ((Na_Oon / Na_Con) (1.0 / 4.0))))
		       (const Na_btfac = (pow ((Na_Ooff / Na_Coff) (1.0 / 4.0))))
		       
		       (const Na_alpha = 150)
		       (const Na_beta  = 3)
		       (const Na_gamma = 150)
		       (const Na_delta = 40)
		       (const Na_epsilon = 1.75)
		       (const Na_zeta = 0.03)
		       (const Na_x1 = 20)
		       (const Na_x2 = -20)
		       (const Na_x3 = 1000000000000.0)
		       (const Na_x4 = -1000000000000.0)
		       (const Na_x5 = 1000000000000.0)
		       (const Na_x6 = -25)))
     ;; end Narsg current

   ))

EOF
))

     (license
      "Copyright 2008-2009 Ivan Raikov and the Okinawa Institute of Science and Technology.

This program 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 3 of the License, or (at
your option) any later version.

This program 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.

A full copy of the GPL license can be found at
<http://www.gnu.org/licenses/>."))))


(define nemo-eggdoc:css (make-parameter #<<EOF
 <!--
      CODE {
            color: #666666;
          }
/*   DT.definition EM { font-weight: bold; font-style: normal; } */

     DT.definition { 
                   background: #eee;
                   color: black;
                   padding: 0.2em 1em 0.2em 0.7em;
                   margin-left: 0.2em;
border: 1px solid #bbc;
                   font-family: "Andale Mono", monospace;
                   /* font-size: 1.2em; */
                   
                 }
     DD {
                   margin-top: 0.8em;
                   margin-bottom: 0.8em;
     }
     DIV.subsection {
                    border-top: 1px solid #448;
                    padding-left: 1em;
                    margin-bottom: 1.2em;
     }
     DIV.subsubsection {
                    border-top: 1px dotted #99c;
                    /* border-left: 1px solid #99c; */
                    padding-left: 1em;
                    margin-bottom: 1.2em;
     }
     DIV.subsubsubsection {
                    border-top: 1px solid #ddf;
                    padding-left: 1em;
                    margin-bottom: 1.2em;
     }

	 DIV.section {
		 margin-bottom: 1.5em;
	 }
	 a:link {
		 color: #336;
	 }
	 a:visited { color: #666; }
	 a:active  { color: #966; }
	 a:hover   { color: #669; }
	 body { margin: 0; padding: 0; background: #fff; color: #000; font: 9pt "Lucida Grande", "Verdana", sans-serif; }
	 H2 {
		 background: #336;
		 color: #fff;
		 padding-top: 0.5em;
		 padding-bottom: 0.5em;
		 padding-left: 16px;
		 margin: 0 0 1em 0;
	}
	UL LI {
		list-style: none;
	}
	TT {
		font-family: "Andale Mono", monospace;
		/* font-size: 1.2em; */
	}
	H3 {
		color: #113;
		margin-bottom: 0.5em;
	}
	H4, H5, H6 {
		color: #113;
		margin-bottom: 1.0em;
	}
	H5 {
		font-weight: normal;
                font-style: italic;
		font-size: 100%;
		margin-top: 1.2em;
	}
	H6 {
		font-weight: bold;
		font-size: 85%;
		margin-top: 1.2em;
	}
     DIV#eggheader {
         text-align: center;
		 float: right;
		 margin-right: 2em;
     }
     DIV#header IMG {
            /* display: block; margin-left: auto; margin-right: auto;  */
            /* float: right; */
            border: none;  /* firefox */
     }
     DIV#footer {
                background: #bbd;
                padding: 0.7em ;
                border-top: 1px solid #cce;
     }
     DIV#footer hr {
                display: none;
     }
     DIV#footer a {
                float: left;
     }
     DIV#revision-history {
         float: right;
     }
     
     DIV#body {
		 margin: 1em 1em 1em 16px;
	 }

     DIV#examples PRE {
       background: #eef;
       padding: 0.1em;
       border: 1px solid #aac;
     }
     PRE#license, DIV#examples PRE {
       padding: 0.5em;
     }
     DIV#examples PRE {
       /* font-size: 85%; */
     }
     PRE { font-family: "Andale Mono", monospace; }
     TABLE {
       background: #eef;
       padding: 0.2em;
       border: 1px solid #aac;
       border-collapse: collapse;
       width: 100%;
     }
     TABLE.symbol-table TD.symbol {
          font-family: "Andale Mono", monospace;
	  vertical-align: top;
          /* font-size: 1.2em; */
     }
     TABLE.symbol-table TD {
          font-family: "Andale Mono", monospace;
	  vertical-align: top;
          /* font-size: 1.2em; */
     }
     TH {
       text-align: left;
       border-bottom: 1px solid #aac;
       padding: 0.25em 0.5em 0.25em 0.5em;
     } 
     TD { padding: 0.25em 0.5em 0.25em 0.5em; }
     -->
EOF
))



(if (eggdoc->html
     doc 
     `( (eggdoc-style . ,(lambda (tag)
			   `("<style type=\"text/css\">" ,(nemo-eggdoc:css)
			     "</style>")))
	
	(comment
	 *macro*
	 . ,(lambda (tag . content) `(p . ,content))) 
			 
	(ebnf-def	; Term definition
	 *macro*
	 . ,(lambda (tag term-def . others)
	      `(dl (dt (b (tt ,term-def)) (tt " ::= "))
		   (dd ,@others))))
			 
	(ebnf-kw	; Keyword
	 *macro*
	 . ,(lambda (tag kw)
	      `((b (tt ,kw)))))
			 
	(ebnf-var	; Variable
	 *macro*
	 . ,(lambda (tag varname) `("{" (i ,varname) "}")))
			 
	(ebnf-choice	; Choice of terms
	 *macro*
	 . ,(lambda (tag . terms)
	      (intersperse terms `(tt " | "))))

	(ebnf-opt	; Optional terms
	 . ,(lambda (tag . terms)
	      `(" [" ,@terms "] ")))
			 
	(ebnf-form	; ( ... ) Form
	 *macro*
	 . ,(lambda (tag . terms)
	      `("( " ,@(intersperse terms " ") " )")))

	(ebnf-compound	; Compound form
	 *macro*
	 . ,(lambda (tag term-def . others)
	      `(dl (dt (b (tt ,term-def)) )
		   (dd ,@others))))
			 
	,@(eggdoc:make-stylesheet doc) ))
    
    (void))