Changeset 30959 in project


Ignore:
Timestamp:
06/04/14 10:24:23 (6 years ago)
Author:
Ivan Raikov
Message:

flsim/signal-diagram: adaptive integration refactoring

Location:
release/4
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • release/4/flsim/trunk/flsim.scm

    r30958 r30959  
    618618exception ConvergenceError
    619619
    620 fun integral (f,x,ys,h,i) =
    621   let open Real
    622       val xn = x+h
    623       val stepper = make_stepper f
    624       fun recur (hs,x,ys,it) =
    625           let
    626               val hf = xn-x
    627           in
    628               if Real.>= (hs, hf)
    629               then (if Real.> (hf,0.0)
    630                     then (let val (stn, etn) = (stepper hf) (x,ys)
    631                           in {tstep=hf, stn=stn} end)
    632                     else {tstep=hs, stn=ys})
    633               else (let val (stn, etn) = (stepper hs) (x,ys)
    634                     in case predictor tol (hs,etn) of
    635                        Left bad => (if Int.<(it,maxit)
    636                                     then recur (bad,x,ys,Int.+(it,1))
    637                                     else raise ConvergenceError)
    638                      | Right good => recur (good,x+hs,stn,0)
    639                    end)
    640           end
    641   in
    642       recur (h,x,ys,0)
    643   end
     620
     621fun secant tol f fg0 guess1 guess0 =
     622    let open Real
     623        val fg1 = f guess1
     624        val newGuess = guess1 - fg1 * (guess1 - guess0) / (fg1 - fg0)
     625        val err =  abs (newGuess - guess1)
     626    in
     627        if (err < tol)
     628        then newGuess
     629        else secant tol f fg1 newGuess guess1
     630    end
     631
     632
     633datatype 'a result = Next of 'a | Root of 'a
     634
     635
     636fun esolver (stepper,evtest) (t,st,tstep) =
     637    let open Real
     638        val (st',e,finterp) = stepper tstep (t,st)
     639    in
     640        case predictor tol (tstep,e) of
     641            Right tstep' =>
     642            if (evtest (st') >= 0.0)
     643            then (let
     644                     val theta   = secant tol (evtest o finterp) (evtest st) 1.0 0.0
     645                     val st''    = finterp (theta+tol)
     646                 in
     647                     Root (t+(theta+tol)*tstep,st'',tstep')
     648                 end)
     649            else Next (t+tstep,st',tstep')
     650          | Left tstep'  =>
     651            esolver (stepper,evtest) (t,st,tstep')
     652    end
     653
     654
     655fun eintegral (f,x,ys,evtest,h,i) =
     656    case esolver (stepper,evtest) (t,st,tstep) of
     657        Next (tn,stn,tstep') =>
     658        ({tn=tn,tstep=tstep',stn=stn})
     659      | Root (tn,stn,tstep') =>
     660        ({tn=tn,stn=stn,tstep=tstep'))
     661
     662fun solver stepper (t,st,tstep) =
     663    let open Real
     664        val (st',e,finterp) = stepper tstep (t,st)
     665    in
     666        case predictor tol (tstep,e) of
     667            Right tstep' =>
     668            (t+tstep,st',tstep')
     669          | Left tstep'  =>
     670            solver (stepper) (t,st,tstep')
     671    end
     672
     673fun integral (f,x,ys,evtest,h,i) =
     674    let
     675        val (tn,stn,tstep') = solver stepper (t,st,tstep)
     676    in
     677        {tn=tn,tstep=tstep',stn=stn}
     678    end
     679
    644680
    645681
  • release/4/signal-diagram/trunk/examples/Izhikevich03.scm

    r30957 r30959  
    2424                 )
    2525
    26             (TRANSIENT subthreshold-eq refractory-eq threshold-detect)
     26            (TRANSIENT subthreshold-eq refractory-eq 'spike threshold-detect)
    2727
    2828            ))
  • release/4/signal-diagram/trunk/signal-diagram.scm

    r30958 r30959  
    255255          )
    256256
     257     (print "fd = " fd)
     258
    257259     (make-sfarrow
    258260      ;; dataflow equations
     
    263265       (lambda (s) (if (prim? fd) outputs (list name)))
    264266       ;; in
    265        (lambda (s) (if (function? fd)
     267       (lambda (s) (begin
     268                     (print "s = " s)
     269                     (print "fd formals = " (function-formals fd))
     270                     (if (function? fd)
    266271                       (lset-intersection eq? (function-formals fd) s)
    267                        s))
     272                       s)))
    268273       ;; out
    269274       (lambda (s) (if (prim? fd) outputs (list name))))
     
    775780     
    776781      ;; in
    777       (lambda (s) (fe-in s))
     782      (lambda (s) (print "sf-actuate in: s = " s) (fe-in s))
    778783     
    779784      ;; out
     
    11321137         (ge-kill (compose (dfe-gen ge) ge-in))
    11331138
    1134          (ee-in   (dfe-in e))
     1139         (ee-in   (dfe-in ee))
    11351140         (ee-out  (compose (dfe-out ee) ee-in))
    11361141         (ee-gen  (compose (dfe-gen ee) ee-in))
     
    11771182
    11781183              (eenv      (map (lambda (s) (cons s s)) (ee-in s)))
    1179               (ecodegen ((sfarrow-codegen f) (ee-in s) eenv ee))
     1184              (ecodegen ((sfarrow-codegen f) (ee-in s) #f eenv ee))
    11801185
    11811186              (fenv      (map (lambda (s) (cons s s)) (fe-in s)))
    1182               (fcodegen ((sfarrow-codegen f) (list e ecompute) (fe-in s) fenv fe))
     1187              (fcodegen ((sfarrow-codegen f) (fe-in s) (list e ecompute (ee-in s)) fenv fe))
    11831188
    11841189              (genv      (map (lambda (s) (cons s s)) (ge-in s)))
    1185               (gcodegen ((sfarrow-codegen g) ev (ge-in s) genv ge))
    1186 
     1190              (gcodegen ((sfarrow-codegen g) (ge-in s) ev genv ge))
    11871191
    11881192              )
     
    12421246
    12431247     ;; signature
    1244      `(TRANSIENT ,(sfarrow-sig f) ,(sfarrow-sig g) ,e)
     1248     `(TRANSIENT ,(sfarrow-sig f) ,(sfarrow-sig g) ,e ,(sfarrow-sig ef))
    12451249
    12461250     ;; children
    1247      `(TRANSIENT ,f ,g)
     1251     `(TRANSIENT ,f ,g ,ef)
    12481252
    12491253     ;; relations
     
    13031307             (rv2 (gensym 'integral))
    13041308             (dfn (gensym 'dfn))
    1305              (evtest (and ev (gensym 'evtest)))
    13061309             )
     1310
    13071311             
    13081312         (lambda (s ev env dfe)
    13091313
    1310            (let* (
     1314           (let* ((evtest (and ev (gensym 'evtest)))
     1315
    13111316                  (idxv      (V:C idx))
    13121317
     
    13381343                                    ))))
    13391344
    1340                (if varh
    1341 
    1342                    ((lambda (l)
    1343                       (if ev
    1344                           (cons
     1345               (if (and varh ev)
     1346                   
     1347                   (let ((ev-select
     1348                          (lambda (x)
     1349                            (let ((yi (find (lambda (y) (equal? x y)) ys)))
     1350                              (if yi
     1351                                  (V:Sub (list-ref yis yi) (V:Var 'yvec))
     1352                                  (select-signal 'ev-select x env))))))
     1353                     
     1354                         (list
    13451355                           (B:Val evtest
    1346                                   (V:Fn `(,x yvec)
    1347                                         (E:Ret (V:Op (cadr ev) (map (lambda (y i) (B:Val y (V:Sub i (V:Var 'yvec)))) ys yis)))
    1348                                         )) l)
    1349                           l))
    1350 
    1351                     (list
    1352                      (B:Val rv1
    1353                             (V:Op 'eintegral
    1354                                   (list (V:Var dfn)
    1355                                         (select-signal 'eintegral1 x env)
    1356                                         (V:Vec (map (lambda (y) (select-signal 'eintegral2 y env)) ys))
    1357                                         (if ev
    1358                                             (V:Op 'SOME (list (V:Var evtest)))
    1359                                             (V:Op 'NONE '()))
    1360                                         tstep
    1361                                         idxv
    1362                                         )))
    1363                      ))
    1364 
     1356                                  (V:Fn `(yvec)
     1357                                        (E:Ret (V:Op (cadr ev) (list (V:Rec (map (lambda (s) `(,s ,(ev-select s))) (caddr ev))))))
     1358                                        ))
     1359
     1360                           (B:Val rv1
     1361                                  (V:Op 'eintegral
     1362                                        (list (V:Var dfn)
     1363                                              (select-signal 'eintegral1 x env)
     1364                                              (V:Vec (map (lambda (y) (select-signal 'eintegral2 y env)) ys))
     1365                                              (V:Var evtest)
     1366                                              tstep
     1367                                              idxv
     1368                                              )))
     1369                           ))
     1370                   
    13651371                   (list
    13661372                    (B:Val rv1
     
    14231429           (RTRANSITION (f g ef eg s)  (sf-rtransition (recur f) (recur g) ef eg s))
    14241430           (TRANSITION (f g ef s)      (sf-transition (recur f) (recur g) ef s))
    1425            (TRANSIENT (f g e)          (sf-transient (recur f) (recur g) (recur e)))
     1431           (TRANSIENT (f g e ef)       (sf-transient (recur f) (recur g) e (recur ef)))
    14261432           (ON (f e)                   (sf-on (recur f) e))
    14271433           (INTEGRAL (x ys h fs)       (sf-integral x ys h fs))
Note: See TracChangeset for help on using the changeset viewer.