Changeset 25981 in project


Ignore:
Timestamp:
02/27/12 10:10:44 (9 years ago)
Author:
Ivan Raikov
Message:

nemo: additions to the nest backend

File:
1 edited

Legend:

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

    r25870 r25981  
    428428      )
    429429
    430   ;; TODO parameters
     430  ;; TODO parameters, time
    431431           
    432432  (for-each (lambda (def)
     
    485485)
    486486
    487 #|
    488 (define (output-init sysname globals state-index-map steady-state-index-map
     487
     488(define (output-parameters sysname const-defs indent indent+)
     489
     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
     550EOF
     551sysname sysname sysname N N N)))
     552
     553
     554
     555(define (output-init sysname state-index-map steady-state-index-map
    489556                     const-defs asgn-eq-defs init-eq-defs rate-eq-defs
    490557                     reaction-eq-defs i-eqs pool-ions perm-ions indent indent+)
    491558
    492     (pp indent ,nl (function y0 = ,(s+ sysname "_init") (Vinit)))
     559 
     560  (let* ((N (length state-index-map))
     561         (init-eqs
     562          (append
     563           
     564           asgn-eq-defs
     565           init-eq-defs
     566           
     567           (map (lambda (pool-ion)
     568                  (let ((n (third pool-ion))
     569                        (b (first pool-ion)))
     570                    (list n b)))
     571                pool-ions)))
     572         
     573         (init-dag
     574          (map (lambda (def)
     575                 (cons (first def) (enum-freevars (second def) '() '())))
     576               init-eqs))
     577         (init-order
     578          (reverse
     579           (topological-sort init-dag
     580                             (lambda (x y) (string=? (->string x) (->string y)))))))
    493581   
    494     (if (not (null? globals)) (pp indent+ (global ,(slp " " globals))))
     582    (pp indent ,nl (,(s+ sysname "::State_::State_") (const Parameters_& p) ": r_(0)"))
    495583   
    496     (pp indent+ ,(expr->string/C++ `(zeros ,(length state-index-map) 1) 'y0))
     584    (pp indent  #\{)
    497585   
    498     (if (member 'v globals)
    499         (let ((vi (lookup-def 'v state-index-map)))
    500           (pp indent+ ,(expr->string/C++ `Vinit 'v))
    501           (if vi (pp indent+ ,(expr->string/C++ 'v (s+ "y0(" vi ")") )))))
    502 
    503 
    504     (let* ((init-eqs
    505             (append
    506              
    507              const-defs
    508              asgn-eq-defs
    509              init-eq-defs
    510              
    511              (map (lambda (pool-ion)
    512                     (let ((n (third pool-ion))
    513                           (b (first pool-ion)))
    514                       (list n b)))
    515                   pool-ions)))
    516 
    517            (init-dag
    518             (map (lambda (def)
    519                    (cons (first def) (enum-freevars (second def) '() '())))
    520                  init-eqs))
    521            (init-order
    522             (reverse
    523              (topological-sort init-dag
    524                (lambda (x y) (string=? (->string x) (->string y)))))))
    525 
    526       (for-each (lambda (n)
    527                   (let ((b  (lookup-def n init-eqs)))
    528                     (if b (pp indent+ ,(expr->string/C++ b (nest-name n))))))
    529                 init-order))
    530 
     586    (pp indent+ ,(expr->string/C++ `(zeros ,N 1) 'y_))
     587
     588    (for-each (lambda (n)
     589                (let ((b  (lookup-def n init-eqs)))
     590                  (if b (pp indent+ ,(expr->string/C++ b (nest-name n))))))
     591              init-order)
     592   
    531593    (for-each (lambda (def)
    532594                (let* ((n  (first def))
    533595                       (ni (lookup-def n state-index-map)))
    534                   (if ni (pp indent+ ,(expr->string/C++ n  (s+ "y0(" ni ")"))))))
     596                  (if ni (pp indent+ ,(expr->string/C++ n  (s+ "y_[" ni "]"))))))
    535597              init-eq-defs)
    536 
     598   
    537599    (if (not (null? steady-state-index-map))
    538600        (begin
    539           (pp indent+ ,(expr->string/C++ `(zeros ,(length steady-state-index-map) 1) 'xs)
    540               ,(expr->string/C++
    541                 `(fsolve ,(s+ #\' sysname "_steadystate" #\' ) xs)
    542                 "[ys,fval,info]"))
    543          
    544          
     601          (gsl-fminimizer sysname N indent+)
    545602          (for-each
    546603           (lambda (def)
     
    549606                    (ni (lookup-def n state-index-map)))
    550607               (if si (begin
    551                         (pp indent+ ,(expr->string/C++ (s+ "ys(" si ")") n))
    552                         (pp indent+ ,(expr->string/C++ (s+ "ys(" si ")") (s+ "y0(" ni ")")))))))
     608                        (pp indent+ ,(expr->string/C++ (s+ "gsl_vector_get(ys, " si ")") n))
     609                        (pp indent+ ,(expr->string/C++ (s+ "gsl_vector_get(ys," si ")")
     610                                                       (s+ "y_[" ni "]")))))))
    553611           rate-eq-defs)
    554 
    555 
     612      (pp indent+ #<<EOF
     613  gsl_vector_free (ys);
     614EOF
     615)
    556616          ))
    557617
    558     (for-each (lambda (def)
    559                 (let ((n (first def)) (b (second def)))
    560                   (if (not (lookup-def n init-eq-defs))
    561                       (pp indent+ ,(expr->string/C++ b n)))))
    562             reaction-eq-defs)
    563 
    564     (for-each (lambda (def)
    565                 (pp indent+ ,(expr->string/C++ (second def) (first def))))
    566               i-eqs)
     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)
     629     
     630      (pp indent  #\})
     631
     632      (pp indent ,nl (,(s+ sysname "::State_::State_") (const Parameters_& p) ": r_(s.r_)"))
    567633   
    568     (pp indent ,nl (end))
    569     (pp indent ,nl (}))
    570 
    571     )
     634      (pp indent  #\{)
     635
     636    (pp indent+ (sprintf #<<EOF
     637    for ( int i = 0 ; i < ~A ; ++i )
     638      y_[i] = s.y_[i];
     639EOF
     640N))
     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); 
     649     for ( size_t i = 0 ; i < ~A ; ++i )
     650       y_[i] = s.y_[i];
     651     r_ = s.r_;
     652     return *this;
     653EOF
     654N))
     655
     656     (pp indent  #\})
     657
     658))
     659
     660
     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)") ))
     666  (pp indent  #\{)
     667
     668  (for-each
     669   (lambda (def)
     670     (let ((n (first def)))
     671       (pp indent+ (sprintf "def<double_t>(d, names::~A,      ~A);" n n))))
     672   const-defs)
     673
     674  (pp indent  #\})
     675
     676  (pp indent ,nl (,(s+ sysname "::Parameters_::set (DictionaryDatum &d))") ))
     677  (pp indent  #\{)
     678
     679  (for-each
     680   (lambda (def)
     681     (let ((n (first def)))
     682       (pp indent+ (sprintf "updateValue<double_t>(d, names::~A,      ~A);" n n))))
     683   const-defs)
     684
     685  (pp indent  #\})
     686
     687
     688#| 
     689  void nest::hh_cond_exp_traub::State_::get(DictionaryDatum &d) const
     690  {
     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]);
     695  }
     696 
     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  }
    572707|#
     708)
    573709
    574710(define (nemo:nest-translator sys . rest)
     
    815951                 ))
    816952
    817 
    818 #|
     953             ;; system parameters
     954             (with-output-to-port output
     955               (lambda ()
     956                 (output-parameters sysname const-defs
     957                                    indent indent+)
     958                 ))
    819959
    820960             ;; initial values function
    821961             (with-output-to-port output
    822962               (lambda ()
    823                  (output-init sysname globals state-index-map steady-state-index-map
     963                 (output-init sysname state-index-map steady-state-index-map
    824964                              const-defs asgn-eq-defs init-eq-defs rate-eq-defs
    825965                              reaction-eq-defs i-eqs pool-ions perm-ions indent indent+)
    826966                 (pp indent ,nl)
    827967                 ))
    828                
    829 |#
    830 
    831            ;; user-defined functions
    832            (let ((define-fn  (make-define-fn)))
    833              (for-each (lambda (fndef)
    834                          (if (not (member (car fndef) builtin-fns))
    835                              (with-output-to-port output
    836                                (lambda ()
    837                                  (apply define-fn (cons indent fndef))
    838                                  (pp indent ,nl)))
    839                              ))
    840                        defuns))
    841          
    842            (close-output-port output))
    843 
     968
     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))
     979             
     980             (close-output-port output))
     981           
    844982           )))
    845983  ))
Note: See TracChangeset for help on using the changeset viewer.