Changeset 26013 in project


Ignore:
Timestamp:
03/01/12 04:01:11 (9 years ago)
Author:
Ivan Raikov
Message:

nemo: further work on the nest backend

Location:
release/4/nemo
Files:
9 edited

Legend:

Unmodified
Added
Removed
  • release/4/nemo/tags/5.0/nemo-nest.scm

    r25858 r26013  
    428428      )
    429429
    430   ;; TODO parameters
     430  ;; TODO parameters, time
    431431           
    432432  (for-each (lambda (def)
     
    815815                 ))
    816816
    817 
    818817#|
    819818
  • release/4/nemo/tags/5.0/nemo-pyparams.scm

    r25858 r26013  
    470470               (lambda (x y) (string=? (->string x) (->string y)))))))
    471471
     472      (pp indent ,nl (Parameters = collections.namedtuple ("'Parameters'" "," []) ))
     473
    472474      (for-each (lambda (x)
    473475                  (pp indent ,(expr->string/python (cadr x) (python-name (car x)))))
  • release/4/nemo/tags/5.1/nemo-nest.scm

    r25871 r26013  
    428428      )
    429429
    430   ;; TODO parameters
     430  ;; TODO parameters, time
    431431           
    432432  (for-each (lambda (def)
     
    815815                 ))
    816816
    817 
    818817#|
    819818
  • release/4/nemo/tags/5.1/nemo-pyparams.scm

    r25871 r26013  
    470470               (lambda (x y) (string=? (->string x) (->string y)))))))
    471471
     472      (pp indent ,nl (Parameters = collections.namedtuple ("'Parameters'" "," []) ))
     473
    472474      (for-each (lambda (x)
    473475                  (pp indent ,(expr->string/python (cadr x) (python-name (car x)))))
  • release/4/nemo/trunk/nemo-core.scm

    r25987 r26013  
    2525
    2626 (make-nemo-core nemo:error nemo:warning
    27                 nemo:env-copy nemo-intern nemo-scoped nemo:quantity?
    28                 nemo:rhs? nemo:conseq? nemo:subst-term nemo:binding? nemo:bind
    29                 eval-nemo-system-decls
    30                 CONST ASGN REACTION RATE PRIM)
     27  nemo:env-copy nemo-intern nemo-scoped nemo:quantity?
     28  nemo:rhs? nemo:conseq? nemo:subst-term nemo:binding? nemo:bind
     29  eval-nemo-system-decls
     30  CONST ASGN REACTION RATE PRIM)
    3131
    3232 (import scheme chicken data-structures ports lolevel extras
  • release/4/nemo/trunk/nemo-matlab.scm

    r25870 r26013  
    679679  (define (cid x)  (second x))
    680680  (define (cn x)   (first x))
    681   (let-optionals rest ((mode 'multiple) (method 'lsode) (filename #f))
     681  (let-optionals rest ((mode 'multiple) (method 'lsode) (filename #t) (dirname "."))
    682682  (match-let ((($ nemo:quantity 'DISPATCH  dis) (environment-ref sys (nemo-intern 'dispatch))))
    683683    (let ((imports  ((dis 'imports)  sys))
     
    687687             (eval-const  (dis 'eval-const))
    688688             (sysname     (matlab-name ((dis 'sysname) sys)))
    689              (prefix      sysname)
    690              (filename    (or filename (s+ sysname ".m")))
     689             (prefix      (->string sysname))
    691690             (deps*       ((dis 'depgraph*) sys))
    692691             (consts      ((dis 'consts)  sys))
     
    907906
    908907           (let ((output  (case mode
    909                             ((single)  (open-output-file filename))
     908                            ((single)  (open-output-file (make-output-fname dirname prefix ".m" filename)))
    910909                            (else #f))))
    911910
     
    913912             
    914913             ;; derivative function
    915              (let ((output1 (or output (open-output-file (s+ prefix "_dy.m")))))
     914             (let ((output1 (or output (open-output-file (make-output-fname dirname prefix "_dy.m")))))
    916915               (with-output-to-port output1
    917916                 (lambda ()
     
    924923
    925924             ;; steady state solver
    926              (let ((output1 (or output (open-output-file (s+ prefix "_steadystate.m")))))
     925             (let ((output1 (or output (open-output-file (make-output-fname dirname prefix "_steadystate.m")))))
    927926               (with-output-to-port output1
    928927                 (lambda ()
     
    933932
    934933             ;; initial values function
    935              (let ((output1 (or output (open-output-file (s+ prefix "_init.m")))))
     934             (let ((output1 (or output (open-output-file (make-output-fname dirname prefix "_init.m")))))
    936935               (with-output-to-port output1
    937936                 (lambda ()
     
    947946             (for-each (lambda (fndef)
    948947                         (if (not (member (car fndef) builtin-fns))
    949                              (let ((output1 (or output (open-output-file (s+ (matlab-name (car fndef)) ".m")))))
     948                             (let ((output1 (or output (open-output-file (make-output-fname dirname (matlab-name (car fndef)) ".m")))))
    950949                               (with-output-to-port output1
    951950                                 (lambda ()
  • release/4/nemo/trunk/nemo-nest.scm

    r25981 r26013  
    2323        (nemo:nest-translator)
    2424
    25         (import scheme chicken utils data-structures lolevel ports srfi-1 srfi-13)
     25        (import scheme chicken utils data-structures lolevel ports extras srfi-1 srfi-13)
    2626       
    2727        (require-extension lolevel matchable strictly-pretty environments
     
    358358
    359359
     360(define (find-locals defs)
     361  (concatenate
     362   (map (lambda (def)
     363          (match def
     364                 (('let bnds body)
     365                  (let ((bexprs (map second bnds)))
     366                    (concatenate (list (map first bnds)
     367                                       (find-locals bexprs )
     368                                       (find-locals (list body))))))
     369                 (('if c t e)      (append (find-locals (list t)) (find-locals (list e))))
     370                 ((s . rest)       (find-locals rest))
     371                 (else (list))))
     372        defs)))
     373
     374
    360375(define (reaction-power sys n)
    361376  (let ((en (environment-ref sys n)))
     
    390405               (body1 (canonicalize-expr/C++ body0))
    391406               (lbs   (enum-bnds body1 (list))))
    392           (if (not (null? lbs)) (pp indent+ (double ,(slp ", " lbs))))
     407          (pp indent+ (double ,rv ";"))
     408          (if (not (null? lbs)) (pp indent+ (double ,(slp ", " lbs) ";")))
    393409          (pp indent+ ,(expr->string/C++ body1 (nest-name rv)))
    394410          (pp indent+ ,(s+ "return " rv ";"))
     
    402418                   pool-ions mcap i-eqs v-eq indent indent+)
    403419
    404   (pp indent ,nl (int ,(s+ sysname "_dynamics")
     420
     421  (pp indent ,nl (int ,(s+ "extern \"C\" " sysname "_dynamics")
    405422                      (,(slp ", " `("double t"
    406423                                    "const double y[]"
     
    430447  ;; TODO parameters, time
    431448           
    432   (for-each (lambda (def)
    433               (let ((n (first def)) )
    434                 (pp indent+
    435                     ,(s+ "// rate equation " n)
    436                     ,(expr->string/C++
    437                       (s+ "y[" (lookup-def n state-index-map) "]") n))))
    438             rate-eq-defs)
    439 
    440   (let* ((eqs
     449
     450  (let* (
     451         (asgn-eq-defs
     452          (map
     453           (lambda (def) (list (first def) (canonicalize-expr/C++ (second def))))
     454           asgn-eq-defs))
     455
     456         (rate-eq-defs
     457          (map
     458           (lambda (def) (list (first def) (canonicalize-expr/C++ (second def))))
     459           rate-eq-defs))
     460
     461         (reaction-eq-defs
     462          (map
     463           (lambda (def) (list (first def) (canonicalize-expr/C++ (second def))))
     464           reaction-eq-defs))
     465
     466         (i-eqs
     467          (map
     468           (lambda (def) (list (first def) (canonicalize-expr/C++ (second def))))
     469           i-eqs))
     470
     471         (v-eq
     472          (and v-eq
     473               (list (first v-eq) (canonicalize-expr/C++ (second v-eq)))))
     474
     475         (eqs
    441476          (append
    442477           
     
    454489                 (cons (first def) (enum-freevars (second def) '() '())))
    455490               eqs))
     491
    456492         (eq-order
    457493          (reverse
    458494           (topological-sort eq-dag
    459                              (lambda (x y) (string=? (->string x) (->string y)))))))
    460           (for-each (lambda (n)
    461                       (let ((b  (lookup-def n eqs)))
    462                         (if b (pp indent+
    463                                   ,(s+ "// equation for " n)
    464                                   ,(expr->string/C++ b (nest-name n))))))
    465                     eq-order))
    466  
    467   (for-each (lambda (def)
    468               (let ((n (s+ "f[" (lookup-def (first def) state-index-map) "]"))
    469                     (b (second def)))
    470                 (pp indent+ ,(s+ "// state " (first def))
    471                     ,(expr->string/C++ b n))))
    472             rate-eq-defs)
    473  
    474   (for-each (lambda (def)
    475               (pp indent+ ,(expr->string/C++ (second def) (first def))))
     495                             (lambda (x y) (string=? (->string x) (->string y))))))
     496
     497         (eq-locals  (find-locals
     498                      (map second
     499                           (if v-eq (cons v-eq (append i-eqs rate-eq-defs eqs))
     500                               (append i-eqs rate-eq-defs eqs)))))
     501         )
     502
     503    (if (not (null? eq-locals))
     504        (pp indent+ (double ,(slp ", " eq-locals) ";")))
     505   
     506    (for-each (lambda (def)
     507                (let ((n (first def)) )
     508                  (pp indent+
     509                      ,(s+ "// rate equation " n)
     510                      ,(expr->string/C++
     511                        (s+ "y[" (lookup-def n state-index-map) "]") n))))
     512              rate-eq-defs)
     513
     514    (for-each (lambda (n)
     515                (let ((b  (lookup-def n eqs)))
     516                  (if b (pp indent+
     517                            ,(s+ "// equation for " n)
     518                            ,(expr->string/C++ b (nest-name n))))))
     519              eq-order)
     520
     521    (for-each (lambda (def)
     522                (let ((n (s+ "f[" (lookup-def (first def) state-index-map) "]"))
     523                      (b (second def)))
     524                  (pp indent+ ,(s+ "// state " (first def))
     525                      ,(expr->string/C++ b n))))
     526              rate-eq-defs)
     527   
     528    (for-each (lambda (def)
     529                (pp indent+ ,(expr->string/C++ (second def) (first def))))
    476530            i-eqs)
    477 
    478   (if v-eq
    479       (let ((vi (lookup-def 'v state-index-map)))
    480         (if vi
    481             (pp indent+ ,(expr->string/C++ (second v-eq) (s+ "f[" vi "]"))))))
    482  
    483   (pp indent+ ,nl ("return GSL_SUCCESS;"))
    484   (pp indent  #\})
    485 )
     531   
     532    (if v-eq
     533        (let ((vi (lookup-def 'v state-index-map)))
     534          (if vi
     535              (pp indent+ ,(expr->string/C++ (second v-eq) (s+ "f[" vi "]"))))))
     536   
     537    (pp indent+ ,nl ("return GSL_SUCCESS;"))
     538    (pp indent  #\})
     539
     540    ))
     541
    486542
    487543
    488544(define (output-parameters sysname const-defs indent indent+)
    489545
    490   (pp indent ,nl (,(s+ sysname "::Parameters_::Parameters_") () ":"))
    491 
    492   (let ((const-exprs
    493          (slp ",\n"
    494           (map (lambda (def)
    495                  (let* ((n  (first def)) (b (second def)))
    496                    (s+ (nest-name n) "  (" (expr->string/C++ b) ")")))
    497                const-defs) )))
    498     (pp indent+ ,const-exprs)
    499    
    500     (pp indent ("{}"))
    501 
    502     ))
    503 
    504 
    505 (define (gsl-fminimizer sysname N indent)
    506 
    507   (pp indent (sprintf #<<EOF
    508 
    509   int minimizer_status, minimizer_iterc;
    510 
    511   const gsl_multimin_fdfminimizer_type *T;
    512   gsl_multimin_fdfminimizer *S;
    513 
    514   gsl_vector *ys;
    515   gsl_multimin_function_fdf F;
    516 
    517   F.f   = &~A_f;
    518   F.df  = &~A_df;
    519   F.fdf = &~A_fdf;
    520   F.n   = ~A;
    521   F.params = &p;
    522 
    523   ys = gsl_vector_alloc (~A);
    524 
    525   for (int i = 0; i < ~A, i++)
    526   {
    527      gsl_vector_set (ys, i, 0.0);
    528   }
    529 
    530   T = gsl_multimin_fdfminimizer_conjugate_fr;
    531   S = gsl_multimin_fdfminimizer_alloc (T, ~A);
    532 
    533   gsl_multimin_fdfminimizer_set (S, &F, ys, 0.01, 1e-4);
    534 
    535   do
    536     {
    537       minimizer_iterc++;
    538       minimizer_status = gsl_multimin_fdfminimizer_iterate (S);
    539 
    540       if (minimizer_status)
    541         break;
    542 
    543       minimizer_status = gsl_multimin_test_gradient (S->gradient, 1e-3);
    544 
    545     }
    546   while (minimizer_status == GSL_CONTINUE && minimizer_iterc < 100);
    547 
    548   gsl_multimin_fdfminimizer_free (S);
    549 
    550 EOF
    551 sysname sysname sysname N N N)))
    552 
     546  (let* ((parameter-defs
     547          (map
     548           (lambda (def) (list (first def) (canonicalize-expr/C++ (second def))))
     549           const-defs))
     550         (parameter-locals  (find-locals (map second parameter-defs))))
     551
     552    (pp indent ,nl (,(s+ sysname "::Parameters_::Parameters_") () ":"))
     553   
     554    (let ((const-exprs
     555           (slp ",\n"
     556                (map (lambda (def)
     557                       (let* ((n  (first def)) (b (second def)))
     558                         (s+ (nest-name n) "  (" (expr->string/C++ b) ")")))
     559                     const-defs) )))
     560     
     561      (if (not (null? parameter-locals))
     562          (pp indent+ (double ,(slp ", " parameter-locals) ";")))
     563
     564      (pp indent+ ,const-exprs)
     565     
     566      (pp indent ("{}"))
     567     
     568      )))
    553569
    554570
    555571(define (output-init sysname state-index-map steady-state-index-map
    556572                     const-defs asgn-eq-defs init-eq-defs rate-eq-defs
    557                      reaction-eq-defs i-eqs pool-ions perm-ions indent indent+)
     573                     reaction-eq-defs i-eqs v-eq pool-ions perm-ions indent indent+)
    558574
    559575 
    560576  (let* ((N (length state-index-map))
     577
     578         (asgn-eq-defs
     579          (map
     580           (lambda (def) (list (first def) (canonicalize-expr/C++ (second def))))
     581           asgn-eq-defs))
     582         
     583         (init-eq-defs
     584          (map
     585           (lambda (def) (list (first def) (canonicalize-expr/C++ (second def))))
     586           init-eq-defs))
     587
     588         (reaction-eq-defs
     589          (map
     590           (lambda (def) (list (first def) (canonicalize-expr/C++ (second def))))
     591           reaction-eq-defs))
     592
     593         (i-eqs
     594          (map
     595           (lambda (def) (list (first def) (canonicalize-expr/C++ (second def))))
     596           i-eqs))
     597
     598         (v-eq
     599          (and v-eq
     600               (list (first v-eq) (canonicalize-expr/C++ (second v-eq)))))
     601         
    561602         (init-eqs
    562603          (append
     
    570611                    (list n b)))
    571612                pool-ions)))
    572          
     613
    573614         (init-dag
    574615          (map (lambda (def)
    575616                 (cons (first def) (enum-freevars (second def) '() '())))
    576617               init-eqs))
     618
    577619         (init-order
    578620          (reverse
    579621           (topological-sort init-dag
    580                              (lambda (x y) (string=? (->string x) (->string y)))))))
     622                             (lambda (x y) (string=? (->string x) (->string y))))))
     623
     624         (init-locals  (find-locals (map second (append init-eqs i-eqs reaction-eq-defs))))
     625         )
    581626   
    582627    (pp indent ,nl (,(s+ sysname "::State_::State_") (const Parameters_& p) ": r_(0)"))
    583628   
    584629    (pp indent  #\{)
    585    
    586     (pp indent+ ,(expr->string/C++ `(zeros ,N 1) 'y_))
     630
     631    (if (not (null? init-locals))
     632        (pp indent+ (double ,(slp ", " init-locals) ";")))
    587633
    588634    (for-each (lambda (n)
     
    610656                                                       (s+ "y_[" ni "]")))))))
    611657           rate-eq-defs)
    612       (pp indent+ #<<EOF
    613   gsl_vector_free (ys);
    614 EOF
    615 )
    616           ))
    617 
    618       (for-each
    619        (lambda (def)
    620          (let ((n (first def)) (b (second def)))
    621            (if (not (lookup-def n init-eq-defs))
    622                (pp indent+ ,(expr->string/C++ b n)))))
    623        reaction-eq-defs)
    624        
    625       (for-each
    626        (lambda (def)
    627          (pp indent+ ,(expr->string/C++ (second def) (first def))))
    628        i-eqs)
     658      (pp indent+ "gsl_vector_free (ys);")
     659      ))
     660
     661    (for-each
     662     (lambda (def)
     663       (let ((n (first def)) (b (second def)))
     664         (if (not (lookup-def n init-eq-defs))
     665             (pp indent+ ,(expr->string/C++ b n)))))
     666     reaction-eq-defs)
     667   
     668    (for-each
     669     (lambda (def)
     670       (pp indent+ ,(expr->string/C++ (second def) (first def))))
     671     i-eqs)
     672
     673    (if v-eq
     674        (let ((vi (lookup-def 'v state-index-map)))
     675          (if vi
     676              (pp indent+ ,(expr->string/C++ (second v-eq) (s+ "y_[" vi "]"))))))
    629677     
    630       (pp indent  #\})
    631 
    632       (pp indent ,nl (,(s+ sysname "::State_::State_") (const Parameters_& p) ": r_(s.r_)"))
    633    
    634       (pp indent  #\{)
    635 
    636     (pp indent+ (sprintf #<<EOF
    637     for ( int i = 0 ; i < ~A ; ++i )
    638       y_[i] = s.y_[i];
    639 EOF
    640 N))
    641      (pp indent  #\})
    642 
    643      (pp indent ,nl (,(s+ sysname "::State_& " sysname "::State_::operator=(const State_& s)")))
    644 
    645      (pp indent  #\{)
    646 
    647     (pp indent+ (sprintf #<<EOF
    648      assert(this != &s); 
     678    (pp indent  #\})
     679   
     680    (pp indent ,nl (,(s+ sysname "::State_::State_") (const Parameters_& p) ": r_(s.r_)"))
     681    (pp indent  #\{)
     682    (pp indent+ (,(sprintf "for ( int i = 0 ; i < ~A ; ++i ) y_[i] = s.y_[i];" N)))
     683    (pp indent  #\})
     684
     685    (pp indent ,nl (,(s+ sysname "::State_& " sysname "::State_::operator=(const State_& s)")))
     686
     687    (pp indent  #\{)
     688
     689
     690    (pp indent+ (,(sprintf #<<EOF
     691   assert(this != &s); 
    649692     for ( size_t i = 0 ; i < ~A ; ++i )
    650693       y_[i] = s.y_[i];
     
    652695     return *this;
    653696EOF
    654 N))
     697N)))
    655698
    656699     (pp indent  #\})
     
    659702
    660703
    661 (define (output-accessors sysname state-index-map steady-state-index-map
    662                           const-defs asgn-eq-defs init-eq-defs rate-eq-defs
    663                           reaction-eq-defs i-eqs pool-ions perm-ions indent indent+)
    664 
    665   (pp indent ,nl (,(s+ sysname "::Parameters_::get (DictionaryDatum &d) const)") ))
     704
     705(define (gsl-fminimizer sysname N indent)
     706
     707  (pp indent (,(sprintf #<<EOF
     708
     709  int minimizer_status, minimizer_iterc;
     710
     711  const gsl_multimin_fdfminimizer_type *T;
     712  gsl_multimin_fdfminimizer *S;
     713
     714  gsl_vector *ys;
     715  gsl_multimin_function_fdf F;
     716
     717  F.f   = &~A_f;
     718  F.df  = &~A_df;
     719  F.fdf = &~A_fdf;
     720  F.n   = ~A;
     721  F.params = &p;
     722
     723  ys = gsl_vector_alloc (~A);
     724
     725  for (int i = 0; i < ~A, i++)
     726  {
     727     gsl_vector_set (ys, i, 0.0);
     728  }
     729
     730  T = gsl_multimin_fdfminimizer_conjugate_fr;
     731  S = gsl_multimin_fdfminimizer_alloc (T, ~A);
     732
     733  gsl_multimin_fdfminimizer_set (S, &F, ys, 0.01, 1e-4);
     734
     735  do
     736    {
     737      minimizer_iterc++;
     738      minimizer_status = gsl_multimin_fdfminimizer_iterate (S);
     739
     740      if (minimizer_status)
     741        break;
     742
     743      minimizer_status = gsl_multimin_test_gradient (S->gradient, 1e-3);
     744
     745    }
     746  while (minimizer_status == GSL_CONTINUE && minimizer_iterc < 100);
     747
     748  gsl_multimin_fdfminimizer_free (S);
     749
     750EOF
     751sysname sysname sysname sysname N N N)))
     752)
     753
     754
     755
     756(define (output-accessors+modifiers
     757         sysname state-index-map const-defs asgn-eq-defs rate-eq-defs
     758         reaction-eq-defs i-eqs v-eq pool-ions perm-ions indent indent+)
     759
     760  (pp indent ,nl (,(sprintf "void ~A::Parameters_::get (DictionaryDatum &d) const)" sysname) ))
    666761  (pp indent  #\{)
    667762
     
    669764   (lambda (def)
    670765     (let ((n (first def)))
    671        (pp indent+ (sprintf "def<double_t>(d, names::~A,      ~A);" n n))))
     766       (pp indent+ (,(sprintf "def<double_t>(d, names::~A, ~A);" n n)))))
    672767   const-defs)
    673768
    674769  (pp indent  #\})
    675770
    676   (pp indent ,nl (,(s+ sysname "::Parameters_::set (DictionaryDatum &d))") ))
     771  (pp indent ,nl (,(sprintf "void ~A::Parameters_::set (DictionaryDatum &d))" sysname) ))
    677772  (pp indent  #\{)
    678773
     
    680775   (lambda (def)
    681776     (let ((n (first def)))
    682        (pp indent+ (sprintf "updateValue<double_t>(d, names::~A,      ~A);" n n))))
     777       (pp indent+ (,(sprintf "updateValue<double_t>(d, names::~A, ~A);" n n)))))
    683778   const-defs)
    684779
    685780  (pp indent  #\})
    686781
    687 
    688 #| 
    689   void nest::hh_cond_exp_traub::State_::get(DictionaryDatum &d) const
     782  (pp indent ,nl (,(sprintf "void ~A::State_::get (DictionaryDatum &d) const)" sysname) ))
     783  (pp indent  #\{)
     784
     785  (for-each
     786   (lambda (def)
     787     (let ((n (first def)) (i (lookup-def (second def) state-index-map)))
     788       (pp indent+ (,(sprintf "def<double_t>(d, names::~A, y_[~A]);" n i)))))
     789   reaction-eq-defs)
     790
     791  (if v-eq
     792     (let ((n (first v-eq)) (i (lookup-def (first v-eq) state-index-map)))
     793       (pp indent+ (,(sprintf "def<double_t>(d, names::~A, y_[~A]);" n i)))))
     794
     795  (pp indent  #\})
     796
     797  (pp indent ,nl (,(sprintf "void ~A::State_::set (const DictionaryDatum &d, const Parameters_&))" sysname) ))
     798  (pp indent  #\{)
     799
     800  (for-each
     801   (lambda (def)
     802     (let ((n (first def)) (i (lookup-def (second def) state-index-map)))
     803       (pp indent+ (,(sprintf "updateValue<double_t>(d, names::~A, y_[~A]);" n i)))))
     804   reaction-eq-defs)
     805
     806  (if v-eq
     807     (let ((n (first v-eq)) (i (lookup-def (first v-eq) state-index-map)))
     808       (pp indent+ (,(sprintf "updateValue<double_t>(d, names::~A, y_[~A]);" n i)))))
     809
     810  (pp indent  #\})
     811
     812)
     813
     814(define (output-buffers+nodes
     815         sysname state-index-map steady-state-index-map
     816         const-defs asgn-eq-defs init-eq-defs rate-eq-defs
     817         reaction-eq-defs i-eqs pool-ions perm-ions indent indent+)
     818
     819  (let ((N (length state-index-map)))
     820
     821  (pp indent ,nl  (,(sprintf #<<EOF
     822~A::Buffers_::Buffers_(~A& n)
     823    : logger_(n),
     824      s_(0),
     825      c_(0),
     826      e_(0)
     827{
     828    // Initialization of the remaining members is deferred to
     829    // init_buffers_().
     830}
     831EOF
     832  sysname sysname)))
     833
     834  (pp indent  ,nl  (,(sprintf #<<EOF
     835~A::Buffers_::Buffers_(const Buffers_&, ~A& n)
     836    : logger_(n),
     837      s_(0),
     838      c_(0),
     839      e_(0)
     840{
     841    // Initialization of the remaining members is deferred to
     842    // init_buffers_().
     843}
     844EOF
     845  sysname sysname)))
     846
     847
     848  (pp indent  ,nl  (,(sprintf #<<EOF
     849~A::~A()
     850    : Archiving_Node(),
     851      P_(),
     852      S_(P_),
     853      B_(*this)
     854{
     855    recordablesMap_.create();
     856}
     857EOF
     858  sysname sysname)))
     859
     860
     861  (pp indent  ,nl  (,(sprintf #<<EOF
     862~A::~A(const ~A& n)
     863    : Archiving_Node(n),
     864      P_(n.P_),
     865      S_(n.S_),
     866      B_(n.B_, *this)
     867{
     868}
     869EOF
     870  sysname sysname sysname)))
     871
     872  (pp indent (,(s+ sysname "::~" sysname "()")))
     873  (pp indent #<<EOF
     874{
     875    // GSL structs only allocated by init_nodes_(), so we need to protect destruction
     876    if ( B_.s_ ) gsl_odeiv_step_free(B_.s_);
     877    if ( B_.c_ ) gsl_odeiv_control_free(B_.c_);
     878    if ( B_.e_ ) gsl_odeiv_evolve_free(B_.e_);
     879}
     880EOF
     881  )
     882
     883  (pp indent ,nl (,(sprintf #<<EOF
     884  void ~A::init_node_(const Node& proto)
     885{
     886    const ~A& pr = downcast<~A>(proto);
     887    P_ = pr.P_;
     888    S_ = pr.S_;
     889}
     890EOF
     891  sysname sysname sysname)))
     892
     893  (pp indent  ,nl (,(sprintf #<<EOF
     894void ~A::init_state_(const Node& proto)
     895{
     896    const ~A& pr = downcast<~A>(proto);
     897    S_ = pr.S_;
     898}
     899EOF
     900  sysname sysname sysname)))
     901
     902  (pp indent ,nl (,(sprintf #<<EOF
     903void ~A::init_buffers_()
     904{
     905    B_.spike_exc_.clear();         
     906    B_.spike_inh_.clear();         
     907    B_.currents_.clear();           
     908    Archiving_Node::clear_history();
     909
     910    B_.logger_.reset();
     911
     912    B_.step_ = Time::get_resolution().get_ms();
     913    B_.IntegrationStep_ = B_.step_;
     914
     915    B_.I_stim_ = 0.0;
     916
     917    static const gsl_odeiv_step_type* T1 = gsl_odeiv_step_rkf45;
     918 
     919    if ( B_.s_ == 0 )
     920      B_.s_ = gsl_odeiv_step_alloc (T1, ~A);
     921    else
     922      gsl_odeiv_step_reset(B_.s_);
     923   
     924    if ( B_.c_ == 0 ) 
     925      B_.c_ = gsl_odeiv_control_y_new (1e-3, 0.0);
     926    else
     927      gsl_odeiv_control_init(B_.c_, 1e-3, 0.0, 1.0, 0.0);
     928   
     929    if ( B_.e_ == 0 ) 
     930      B_.e_ = gsl_odeiv_evolve_alloc(~A);
     931    else
     932      gsl_odeiv_evolve_reset(B_.e_);
     933 
     934    B_.sys_.function  = ~A_dynamics;
     935    B_.sys_.jacobian  = 0;
     936    B_.sys_.dimension = ~A;
     937    B_.sys_.params    = reinterpret_cast<void*>(this);
     938}
     939EOF
     940  sysname N N sysname N)))
     941
     942  (let ((iv (lookup-def 'v state-index-map)))
     943    (if iv
     944        (pp indent ,nl (,(sprintf #<<EOF
     945void ~A::calibrate()
     946{
     947    B_.logger_.init(); 
     948    V_.RefractoryCounts_ = 20;
     949    V_.U_old_ = S_.y_[~A];
     950}
     951EOF
     952  sysname iv)))
     953  ))
     954
     955))
     956
     957
     958(define (output-update
     959         sysname state-index-map steady-state-index-map
     960         const-defs asgn-eq-defs init-eq-defs rate-eq-defs
     961         reaction-eq-defs i-eqs pool-ions perm-ions indent indent+)
     962 
     963  (let ((vi (lookup-def 'v state-index-map)))
     964   
     965    (pp indent  (,(sprintf #<<EOF
     966void ~A::update(Time const & origin, const long_t from, const long_t to)
     967EOF
     968  sysname)))
     969
     970    (pp indent  #\{)
     971
     972    (pp indent+ (,(sprintf #<<EOF
     973assert(to >= 0 && (delay) from < Scheduler::get_min_delay());
     974    assert(from < to);
     975
     976    for ( long_t lag = from ; lag < to ; ++lag )
     977      {
     978   
     979        double tt = 0.0 ; //it's all relative!
     980        V_.U_old_ = S_.y_[~A];
     981
     982   
     983        // adaptive step integration
     984        while (tt < B_.step_)
     985        {
     986          const int status = gsl_odeiv_evolve_apply(B_.e_, B_.c_, B_.s_,
     987                                 &B_.sys_,              // system of ODE
     988                                 &tt,                   // from t...
     989                                  B_.step_,             // ...to t=t+h
     990                                 &B_.IntegrationStep_ , // integration window (written on!)
     991                                  S_.y_);               // neuron state
     992
     993          if ( status != GSL_SUCCESS )
     994            throw GSLSolverFailure(get_name(), status);
     995        }
     996
     997        // sending spikes: crossing 0 mV, pseudo-refractoriness and local maximum...
     998        // refractory?
     999        if (S_.r_)
     1000          {
     1001            --S_.r_;
     1002          }
     1003        else
     1004          {
     1005            // (threshold && maximum)
     1006            if (S_.y_[~A] >= P_.V_T + 30. && V_.U_old_ > S_.y_[~A])
     1007              {
     1008                S_.r_ = V_.RefractoryCounts_;
     1009               
     1010                set_spiketime(Time::step(origin.get_steps()+lag+1));
     1011               
     1012                SpikeEvent se;
     1013                network()->send(*this, se, lag);
     1014              }
     1015          }
     1016   
     1017        // set new input current
     1018        B_.I_stim_ = B_.currents_.get_value(lag);
     1019
     1020        // log state data
     1021        B_.logger_.record_data(origin.get_steps() + lag);
     1022
     1023      }
     1024EOF
     1025  vi vi vi)))
     1026 
     1027  (pp indent  #\})
     1028
     1029))
     1030
     1031
     1032(define (output-handle
     1033         sysname state-index-map steady-state-index-map
     1034         const-defs asgn-eq-defs init-eq-defs rate-eq-defs
     1035         reaction-eq-defs i-eqs pool-ions perm-ions indent indent+)
     1036 
     1037  (pp indent  ,nl (,(sprintf #<<EOF
     1038void ~A::handle(SpikeEvent & e)
    6901039  {
    691     def<double>(d, names::V_m   , y_[V_M]); // Membrane potential
    692     def<double>(d,names::Act_m  , y_[HH_M]);
    693     def<double>(d,names::Act_h  , y_[HH_H]);
    694     def<double>(d,names::Inact_n, y_[HH_N]);
     1040    assert(e.get_delay() > 0);
     1041
     1042    if(e.get_weight() > 0.0)
     1043      {
     1044        B_.spike_exc_.add_value(e.get_rel_delivery_steps(network()->get_slice_origin()),
     1045                                e.get_weight() * e.get_multiplicity());
     1046      }
     1047    else
     1048      {
     1049        // add with negative weight, ie positive value, since we are changing a
     1050        // conductance
     1051        B_.spike_inh_.add_value(e.get_rel_delivery_steps(network()->get_slice_origin()),
     1052                                -e.get_weight() * e.get_multiplicity());
     1053      }
    6951054  }
     1055EOF
     1056  sysname)))
     1057
     1058  (pp indent  ,nl (,(sprintf #<<EOF
     1059void ~A::handle(CurrentEvent& e)
     1060  {
     1061    assert(e.get_delay() > 0);
     1062
     1063    const double_t c=e.get_current();
     1064    const double_t w=e.get_weight();
     1065
     1066    B_.currents_.add_value(e.get_rel_delivery_steps(network()->get_slice_origin()),
     1067                        w *c);
     1068  }
     1069
     1070void ~A::handle(DataLoggingRequest& e)
     1071  {
     1072    B_.logger_.handle(e);
     1073  }
     1074EOF
     1075  sysname sysname)))
     1076
     1077)
     1078
     1079
     1080(define (output-recordables
     1081         sysname state-index-map steady-state-index-map
     1082         const-defs asgn-eq-defs init-eq-defs rate-eq-defs
     1083         reaction-eq-defs i-eqs pool-ions perm-ions indent indent+)
    6961084 
    697   void nest::hh_cond_exp_traub::State_::set(const DictionaryDatum& d, const Parameters_&)
    698   {
    699     updateValue<double>(d, names::V_m   , y_[V_M]);
    700     updateValue<double>(d,names::Act_m  , y_[HH_M]);
    701     updateValue<double>(d,names::Act_h  , y_[HH_H]);
    702     updateValue<double>(d,names::Inact_n, y_[HH_N]);
    703 
    704     if ( y_[HH_M] < 0 || y_[HH_H] < 0 || y_[HH_N] < 0 )
    705         throw BadProperty("All (in)activation variables must be non-negative.");
    706   }
    707 |#
    708 )
     1085  (pp indent ,nl (,(sprintf "RecordablesMap<~A> ~A::recordablesMap_;" sysname sysname)))
     1086  (pp indent (,(sprintf "template <> void RecordablesMap<~A>::create()" sysname)))
     1087  (pp indent #\{)
     1088  (for-each
     1089   (lambda (def)
     1090     (let ((n (first def)) (i (second def)))
     1091       (pp indent+ (,(sprintf "insert_(names::~A, &~A::get_y_elem_<~A>);"
     1092                              n sysname i)))))
     1093   state-index-map)
     1094  (pp indent #\})
     1095  )
     1096
     1097
     1098(define (output-prelude sysname indent)
     1099  (pp indent (,(sprintf "#include \"~A.h\"
     1100#include \"exceptions.h\"
     1101#include \"network.h\"
     1102#include \"dict.h\"
     1103#include \"integerdatum.h\"
     1104#include \"doubledatum.h\"
     1105#include \"dictutils.h\"
     1106#include \"numerics.h\"
     1107#include <limits>
     1108
     1109#include \"universal_data_logger_impl.h\"
     1110
     1111#include <iomanip>
     1112#include <iostream>
     1113#include <cstdio>
     1114
     1115#include <gsl/gsl_errno.h>
     1116#include <gsl/gsl_matrix.h>
     1117#include <gsl/gsl_sf_exp.h>
     1118"  sysname))))
    7091119
    7101120(define (nemo:nest-translator sys . rest)
     
    7131123  (define (cn x)   (first x))
    7141124
    715   (let-optionals rest ((filename #f))
     1125  (let-optionals rest ((dirname "."))
    7161126  (match-let ((($ nemo:quantity 'DISPATCH  dis) (environment-ref sys (nemo-intern 'dispatch))))
    7171127    (let ((imports  ((dis 'imports)  sys))
     
    7211131             (eval-const  (dis 'eval-const))
    7221132             (sysname     (nest-name ((dis 'sysname) sys)))
    723              (prefix      sysname)
    724              (filename    (or filename (s+ sysname ".cpp")))
     1133             (prefix      (->string sysname))
    7251134             (deps*       ((dis 'depgraph*) sys))
    7261135             (consts      ((dis 'consts)  sys))
     
    8981307                                                    (list (+ 1 (first ax))
    8991308                                                          (cons `(,st-name ,(first ax)) (second ax)))))
    900                                                 (list 1 (list))
     1309                                                (list 0 (list))
    9011310                                                (if (member 'v globals)
    9021311                                                    (cons (list 'v) rate-eq-defs)
     
    9111320                                                                     (cons `(,st-name ,(first ax)) (second ax)))
    9121321                                                               ax)))
    913                                                        (list 1 (list))
     1322                                                       (list 0 (list))
    9141323                                                       rate-eq-defs)))
    9151324                                        (second acc)))
     
    9401349            pool-ions)
    9411350
    942            (let ((output  (open-output-file filename)))
     1351           (let ((cpp-output  (open-output-file (make-output-fname dirname prefix ".cpp"))))
     1352
     1353             ;; include files and other prelude definitions
     1354             
     1355             (with-output-to-port cpp-output
     1356               (lambda ()
     1357                 (output-prelude sysname indent)
     1358                 ))
     1359
     1360             ;; user-defined functions
     1361             (let ((define-fn  (make-define-fn)))
     1362               (for-each (lambda (fndef)
     1363                           (if (not (member (car fndef) builtin-fns))
     1364                               (with-output-to-port cpp-output
     1365                                 (lambda ()
     1366                                   (apply define-fn (cons indent fndef))
     1367                                   (pp indent ,nl)))
     1368                               ))
     1369                         defuns))
     1370
     1371             ;;  system recordables
     1372             (with-output-to-port cpp-output
     1373               (lambda ()
     1374                 (output-recordables
     1375                  sysname state-index-map steady-state-index-map
     1376                  const-defs asgn-eq-defs init-eq-defs rate-eq-defs
     1377                  reaction-eq-defs i-eqs pool-ions perm-ions indent indent+)
     1378                 (pp indent ,nl)
     1379                 ))
     1380
    9431381             
    9441382             ;; derivative function
    945              (with-output-to-port output
     1383             (with-output-to-port cpp-output
    9461384               (lambda ()
    9471385                 (output-dy sysname globals state-index-map
     
    9521390
    9531391             ;; system parameters
    954              (with-output-to-port output
     1392             (with-output-to-port cpp-output
    9551393               (lambda ()
    9561394                 (output-parameters sysname const-defs
     
    9591397
    9601398             ;; initial values function
    961              (with-output-to-port output
     1399             (with-output-to-port cpp-output
    9621400               (lambda ()
    9631401                 (output-init sysname state-index-map steady-state-index-map
    9641402                              const-defs asgn-eq-defs init-eq-defs rate-eq-defs
    965                               reaction-eq-defs i-eqs pool-ions perm-ions indent indent+)
     1403                              reaction-eq-defs i-eqs v-eq
     1404                              pool-ions perm-ions indent indent+)
    9661405                 (pp indent ,nl)
    9671406                 ))
    9681407
    969              ;; user-defined functions
    970              (let ((define-fn  (make-define-fn)))
    971                (for-each (lambda (fndef)
    972                            (if (not (member (car fndef) builtin-fns))
    973                                (with-output-to-port output
    974                                  (lambda ()
    975                                    (apply define-fn (cons indent fndef))
    976                                    (pp indent ,nl)))
    977                                ))
    978                          defuns))
     1408             ;; accessors and modifiers for states and parameters
     1409             (with-output-to-port cpp-output
     1410               (lambda ()
     1411                 (output-accessors+modifiers
     1412                  sysname state-index-map
     1413                  const-defs asgn-eq-defs rate-eq-defs
     1414                  reaction-eq-defs i-eqs v-eq pool-ions perm-ions
     1415                  indent indent+)
     1416                 (pp indent ,nl)
     1417                 ))
     1418
    9791419             
    980              (close-output-port output))
     1420             ;; buffer and node initialization
     1421             (with-output-to-port cpp-output
     1422               (lambda ()
     1423                 (output-buffers+nodes
     1424                  sysname state-index-map steady-state-index-map
     1425                  const-defs asgn-eq-defs init-eq-defs rate-eq-defs
     1426                  reaction-eq-defs i-eqs pool-ions perm-ions indent indent+)
     1427                 (pp indent ,nl)
     1428                 ))
     1429
     1430             ;; state update
     1431             (with-output-to-port cpp-output
     1432               (lambda ()
     1433                 (output-update
     1434                  sysname state-index-map steady-state-index-map
     1435                  const-defs asgn-eq-defs init-eq-defs rate-eq-defs
     1436                  reaction-eq-defs i-eqs pool-ions perm-ions indent indent+)
     1437                 (pp indent ,nl)
     1438                 ))
     1439
     1440             ;; spike handler
     1441             (with-output-to-port cpp-output
     1442               (lambda ()
     1443                 (output-handle
     1444                  sysname state-index-map steady-state-index-map
     1445                  const-defs asgn-eq-defs init-eq-defs rate-eq-defs
     1446                  reaction-eq-defs i-eqs pool-ions perm-ions indent indent+)
     1447                 (pp indent ,nl)
     1448                 ))
     1449
     1450             
     1451             (close-output-port cpp-output))
    9811452           
    9821453           )))
  • release/4/nemo/trunk/nemo-utils.scm

    r25870 r26013  
    2222       
    2323 (lookup-def enum-bnds enum-freevars sum
    24              if-convert let-enum let-elim let-lift
    25              s+ sw+ slp nl spaces ppf
    26              transitions-graph state-conseqs
    27              differentiate simplify distribute)
    28 
    29  (import scheme chicken data-structures srfi-1 srfi-13)
     24  if-convert let-enum let-elim let-lift
     25  s+ sw+ slp nl spaces ppf
     26  transitions-graph state-conseqs
     27  differentiate simplify distribute
     28  make-output-fname)
     29
     30 (import scheme chicken data-structures files srfi-1 srfi-13)
    3031
    3132 (require-extension matchable strictly-pretty
     
    258259    (list n conseqs1)))
    259260
     261
     262(define (make-output-fname dirname sysname suffix . rest)
     263  (let-optionals rest ((x #t))
     264    (and x
     265         (if (string? x) x
     266             (let ((fname (s+ sysname suffix)))
     267               (or (and dirname (make-pathname dirname fname)) fname)))
     268         )))
     269
     270
    260271;;    `(+ - * / pow neg abs atan asin acos sin cos exp ln
    261272;;      sqrt tan cosh sinh tanh hypot gamma lgamma log10 log2 log1p ldexp cube
  • release/4/nemo/trunk/nemo.scm

    r25870 r26013  
    2121(import scheme chicken srfi-1)
    2222
    23 (require-extension nemo-core nemo-macros nemo-hh nemo-vclamp)
     23(require-extension nemo-core nemo-macros nemo-hh nemo-vclamp nemo-utils)
    2424(require-extension matchable ssax sxml-transforms sxpath sxpath-lolevel
    2525                   environments getopt-long)
     
    9595          `(
    9696            (nest
    97              "write NEST output to given file (default: <model-name>.cpp)"
    98              (value (optional FILENAME)))
    99                      
     97             "write NEST output files <model-name>.cpp and <model-name>.h in the given directory (default: .)"
     98             (value (optional DIRNAME)))
    10099            )
    101100          `())
     
    111110    ,@(if nemo-matlab?
    112111          `((matlab
    113              "write MATLAB output")
     112             "write MATLAB output in the given directory (default: .)"
     113             (value (optional DIRNAME)))
    114114
    115115            (octave
     
    281281  (if nemo-nest?
    282282      (lambda (options model)
    283         (nemo:nest-translator model (lookup-def 'filename options)))
     283        (nemo:nest-translator model (lookup-def 'dirname options)))
    284284      (lambda (options model)
    285285        (void))))
     
    288288  (if nemo-matlab?
    289289      (lambda (options model)
    290         (nemo:matlab-translator model))
     290        (nemo:matlab-translator model #f (lookup-def 'dirname options)))
    291291      (lambda (options model)
    292292        (void))))
     
    305305        (nemo:octave-translator model
    306306                                (lookup-def 'method options)
    307                                 (lookup-def 'filename options)))
     307                                (lookup-def 'filename options)
     308                                (lookup-def 'dirname options)))
    308309      (lambda (options model)
    309310        (void))))
     
    886887        (for-each
    887888         (lambda (operand model)
    888            (let* (
    889                   (make-fname  (lambda (operand x suffix)
    890                                  (and x (let ((fname (pathname-strip-extension operand)))
    891                                           (if (string? x) x (s+ fname suffix))))))
    892                  
    893                   (sxml-fname   (make-fname operand (opt 'sxml) ".sxml"))
    894                   (xml-fname    (make-fname operand (opt 'xml) ".xml"))
    895                   (mod-fname    (make-fname operand (opt 'nmodl) ".mod"))
    896                   (octave-fname (and (string? (opt 'octave)) (make-fname operand (opt 'octave) ".m") ))
    897                   (nest-fname   (and (string? (opt 'nest)) (make-fname operand (opt 'nest) ".cpp") ))
    898                   (vclamp-ses-fname    (make-fname operand (opt 'vclamp-hoc) ".ses"))
    899                   (vclamp-octave-fname (make-fname operand (opt 'vclamp-octave) "_vclamp.m"))
     889
     890           (match-let ((($ nemo:quantity 'DISPATCH  dis) (environment-ref model (nemo-intern 'dispatch))))
     891
     892            (let* ((sysname      ((dis 'sysname) model))
     893                   (dirname      (pathname-directory operand))
     894                   (sxml-fname   (make-output-fname dirname sysname  ".sxml" (opt 'sxml) ))
     895                   (xml-fname    (make-output-fname dirname sysname ".xml"  (opt 'xml) ))
     896                   (mod-fname    (make-output-fname dirname sysname ".mod"  (opt 'nmodl) ))
     897                   (vclamp-ses-fname    (make-output-fname dirname sysname ".ses" (opt 'vclamp-hoc) ))
     898                   (vclamp-octave-fname (make-output-fname dirname sysname "_vclamp.m" (opt 'vclamp-octave) ))
    900899                 
    901900                  (nest           (opt 'nest))
     
    939938                                                model))))
    940939             
    941              (if octave (model->octave `((filename  . ,octave-fname)
    942                                          (method  . ,octave-method))
     940             (if octave (model->octave `((filename  . ,(or (and (string? octave) (pathname-file octave)) octave))
     941                                         (dirname   . ,(or (and (string? octave) (pathname-directory octave)) dirname))
     942                                         (method    . ,octave-method))
    943943                                       model))
    944944             
    945              (if matlab (model->matlab `() model))
     945             (if matlab (model->matlab `((dirname . ,(or (and (string? matlab) matlab) dirname))) model))
    946946             
    947              (if nest (model->nest `((filename . ,nest-fname)) model))
    948              
     947             (if nest (model->nest `((dirname . ,(or (and (string? nest) nest) dirname))) model))
    949948             
    950949             (if vclamp-hoc (model->vclamp-hoc `((filename . ,vclamp-ses-fname)
     
    954953                                                       )
    955954                                                     model))
     955
    956956             
    957              ))
     957             )))
    958958         operands models)
    959 
    960       (let ((pyparams       (opt 'pyparams)))
    961         (if pyparams
    962             (let ((pyparams-fname
    963                    (or (and (string? (opt 'pyparams))
    964                             (make-fname operand (opt 'pyparams) ".py") )
    965                        "pyparams.py")))
    966               (models->pyparams `((filename . ,pyparams-fname)) models))))
     959       
     960        (let ((pyparams (opt 'pyparams)))
     961          (if pyparams
     962              (let ((pyparams-fname
     963                     (or (and (string? pyparams) pyparams)
     964                         "pyparams.py")))
     965                (models->pyparams `((filename . ,pyparams-fname)) models))))
    967966
    968967      )))
Note: See TracChangeset for help on using the changeset viewer.