Changeset 7961 in project


Ignore:
Timestamp:
01/27/08 12:43:13 (12 years ago)
Author:
Ivan Raikov
Message:

Some changes to the solver interface for compatibility with other Scheme implementations.

Location:
ode/trunk
Files:
5 edited

Legend:

Unmodified
Added
Removed
  • ode/trunk/abm4.scm

    r7924 r7961  
    123123       (for-each
    124124        (lambda (x)
    125           (let-values (((df sym rhs)  (match x ((sym rhs) (values #f sym rhs))
    126                                              ((df sym rhs) (values df sym rhs)))))
    127             (let ((v     (eval-rhs rhs))
    128                   (svec  (solve-env-ref sym)))
    129               (cond ((is-state? sym)  (case df
    130                                         ((q)  (svec-value-set! svec v))
    131                                         ((d)  (svec-deriv-set! svec v))))
    132                     ((is-asgn? sym)   (begin
    133                                         (eval-env-set! sym v)
    134                                         (if (vector? svec)
    135                                             (svec-value-set! svec v)
    136                                             (set-car! svec v))))))))
    137         lst))))
     125          (let ((xrec (if (> (length x) 2) x (cons #f x))))
     126            (let ((df   (first xrec))
     127                  (sym  (second xrec))
     128                  (rhs  (third xrec)))
     129              (let ((v     (eval-rhs rhs))
     130                    (svec  (solve-env-ref sym)))
     131                (cond ((is-state? sym)  (case df
     132                                          ((q)  (svec-value-set! svec v))
     133                                          ((d)  (svec-deriv-set! svec v))))
     134                      ((is-asgn? sym)   (begin
     135                                          (eval-env-set! sym v)
     136                                          (if (vector? svec)
     137                                              (svec-value-set! svec v)
     138                                              (set-car! svec v)))))))))
     139          lst))))
    138140 
    139141
     
    359361
    360362            (begin
    361               (values it t tstep))))))
     363              (list it t tstep))))))
    362364       
    363365  ;; Predictor-corrector
     
    368370    (eval-derivs!)
    369371
    370     (let-values (((it t tstep)  (rk5 it t tstep #t)))
    371 
    372       (let loop ((it it) (t t) (tstep tstep))
    373 
    374         (for-each (lambda (g)
    375                     (if (not (eval-rhs g))
    376                         (numerror "guard failed: " g))) guards)
    377 
    378         (if (t<tstop tstep t tstop)
    379             (begin
    380               (if (> (* tstep (- (+ t tstep) tstop)) 0.0)
    381                   (begin
    382                     (set! tstep (- tstop t))
    383                     (eval-env-set! dt tstep)))
    384 
    385               ;; Adams-Bashforth predictor
    386               ((ode-runtime 'solve-env-state-fold)
    387                (lambda (sym svec i+hstlist)
    388                  (let ((i        (car i+hstlist))
    389                        (hstlist  (cdr i+hstlist)))
    390                    (let* ((sthist   (car hstlist))
    391                           (pri      (cons (svec-deriv svec) (map cdr (stack->list sthist))))
    392                           (predval  (apply ab-predictor
    393                                            (cons (car (stack-peek sthist))
    394                                                  (cons tstep pri)))))
     372    (let ((rk5rec (rk5 it t tstep #t)))
     373      (let ((it     (first rk5rec))
     374            (t      (second rk5rec))
     375            (tstep  (third rk5rec)))
     376
     377        (let loop ((it it) (t t) (tstep tstep))
     378
     379          (for-each (lambda (g)
     380                      (if (not (eval-rhs g))
     381                          (numerror "guard failed: " g))) guards)
     382
     383          (if (t<tstop tstep t tstop)
     384              (begin
     385                (if (> (* tstep (- (+ t tstep) tstop)) 0.0)
     386                    (begin
     387                      (set! tstep (- tstop t))
     388                      (eval-env-set! dt tstep)))
     389               
     390                ;; Adams-Bashforth predictor
     391                ((ode-runtime 'solve-env-state-fold)
     392                 (lambda (sym svec i+hstlist)
     393                   (let ((i        (car i+hstlist))
     394                         (hstlist  (cdr i+hstlist)))
     395                     (let* ((sthist   (car hstlist))
     396                            (pri      (cons (svec-deriv svec) (map cdr (stack->list sthist))))
     397                            (predval  (apply ab-predictor
     398                                             (cons (car (stack-peek sthist))
     399                                                   (cons tstep pri)))))
    395400                       (vector-set! predvect i predval)
    396401                       (svec-value-set! svec predval)
    397402                       (eval-env-set! sym predval)
    398403                       (cons (fx+ 1 i) (cdr hstlist)))))
    399               (cons 0 history))
    400              
    401               (eval-env-set! indep t)
    402               (eval-env-set! dt tstep)
    403               (eval-derivs!)
    404              
    405               ;; Adams-Moulton corrector
    406               ((ode-runtime 'solve-env-state-fold)
    407               (lambda (sym svec i+hstlist)
    408                  (let ((i        (car i+hstlist))
    409                       (hstlist  (cdr i+hstlist)))
    410                    (let* ((sthist   (car hstlist))
    411                           (pri      (cons (svec-deriv svec) (map cdr (stack->list sthist))))
    412                           (val      (apply am-corrector
    413                                            (cons (car (stack-peek sthist))
    414                                                  (cons tstep pri)))))
     404                (cons 0 history))
     405               
     406                (eval-env-set! indep t)
     407                (eval-env-set! dt tstep)
     408                (eval-derivs!)
     409               
     410                ;; Adams-Moulton corrector
     411                ((ode-runtime 'solve-env-state-fold)
     412                (lambda (sym svec i+hstlist)
     413                   (let ((i        (car i+hstlist))
     414                        (hstlist  (cdr i+hstlist)))
     415                     (let* ((sthist   (car hstlist))
     416                            (pri      (cons (svec-deriv svec) (map cdr (stack->list sthist))))
     417                            (val      (apply am-corrector
     418                                             (cons (car (stack-peek sthist))
     419                                                   (cons tstep pri)))))
    415420                       (let ((predval (vector-ref predvect i)))
    416421                         (if (not (= val 0.0))
     
    421426                           (svec-value-set! svec cval))
    422427                         (cons (fx+ 1 i) (cdr hstlist))))))
    423                (cons 0 history))
    424              
    425               (let ((err (maxerror)))
    426                 (let ((hierr (hierror? t tstep err))
    427                       (loerr (lowerror? t tstep err)))
    428                   (cond (hierr
    429                          (begin
    430                            (set! tstep (* tstep (max m (sqrt (/ relmax (abs hierr))))))
    431                            ((ode-runtime 'solve-env-state-fold)
    432                             (lambda (sym svec hstlist)
    433                               (let* ((sthist   (car hstlist))
    434                                      (pri      (cdr (stack-peek sthist)))
    435                                      (val      (car (stack-peek sthist))))
    436                                 (eval-env-set! sym val)
    437                                 (svec-value-set! svec val)
    438                                 (svec-deriv-set! svec pri))
    439                               (cdr hstlist))
    440                             history)
    441                            (let-values (((it t tstep) (rk5 it (eval-env-ref indep) tstep #f)))
    442                                        (loop it t tstep))))
    443                        
    444                         (else
    445                          (begin
    446                            
    447 
    448                            ;; update the derivative and value history
    449                            ((ode-runtime 'solve-env-state-fold)
    450                             (lambda (sym svec hstlist)
    451                               (let ((sthist  (car hstlist)))
    452                                 (let ((sthist1 (stack-push! sthist
    453                                                             (cons (svec-value svec)
    454                                                                   (svec-deriv svec)))))
    455                                   (let ((dep (stack-depth sthist1)))
    456                                     (if (fx<= history-len dep) (stack-cut! sthist1 history-len dep)))
    457                                   (cdr hstlist))))
    458                             history)
    459                            ;; update independent variable
    460                            (solve-env-set! indep svec-value-idx (+ t tstep))
    461                            (printq)
    462 
    463                            (if (and loerr (positive? loerr))
     428                 (cons 0 history))
     429               
     430                (let ((err (maxerror)))
     431                  (let ((hierr (hierror? t tstep err))
     432                        (loerr (lowerror? t tstep err)))
     433                    (cond (hierr
     434                           (begin
     435                             (set! tstep (* tstep (max m (sqrt (/ relmax (abs hierr))))))
     436                             ((ode-runtime 'solve-env-state-fold)
     437                              (lambda (sym svec hstlist)
     438                                (let* ((sthist   (car hstlist))
     439                                       (pri      (cdr (stack-peek sthist)))
     440                                       (val      (car (stack-peek sthist))))
     441                                  (eval-env-set! sym val)
     442                                  (svec-value-set! svec val)
     443                                  (svec-deriv-set! svec pri))
     444                                (cdr hstlist))
     445                              history)
     446                             (apply loop (rk5 it (eval-env-ref indep) tstep #f))))
     447                         
     448                          (else
     449                           (begin
     450                             
     451                             
     452                             ;; update the derivative and value history
     453                             ((ode-runtime 'solve-env-state-fold)
     454                              (lambda (sym svec hstlist)
     455                                (let ((sthist  (car hstlist)))
     456                                  (let ((sthist1 (stack-push! sthist
     457                                                              (cons (svec-value svec)
     458                                                                    (svec-deriv svec)))))
     459                                    (let ((dep (stack-depth sthist1)))
     460                                      (if (fx<= history-len dep) (stack-cut! sthist1 history-len dep)))
     461                                    (cdr hstlist))))
     462                              history)
     463                             ;; update independent variable
     464                             (solve-env-set! indep svec-value-idx (+ t tstep))
     465                             (printq)
     466                             
     467                             (if (and loerr (positive? loerr))
    464468                                 (set! tstep (* tstep S (min M (sqrt (/ relmax (abs loerr)))))))
    465                            
    466                            (eval-env-set! indep t)
    467                            (eval-env-set! dt tstep)
    468                            (eval-derivs!)
    469 
    470                            (loop (fx+ 1 it) (+ t tstep) tstep)))))))
    471             (values it t tstep (list))))))
     469                            
     470                             (eval-env-set! indep t)
     471                             (eval-env-set! dt tstep)
     472                             (eval-derivs!)
     473                             
     474                             (loop (fx+ 1 it) (+ t tstep) tstep))))))))
     475            (list it t tstep (list))))))
    472476 
    473477  (if debug (print "abm4: " (ode-runtime 'debug)))   
  • ode/trunk/euler.scm

    r7924 r7961  
    7878       (fold
    7979        (lambda (x i)
    80           (let-values (((df sym rhs)  (match x ((sym rhs) (values #f sym rhs))
    81                                              ((df sym rhs) (values df sym rhs)))))
     80          (let ((xrec (if (> (length x) 2) x (cons #f x))))
     81            (let ((df   (first xrec))
     82                  (sym  (second xrec))
     83                  (rhs  (third xrec)))
    8284            (let ((v     (eval-rhs rhs))
    8385                  (svec  (solve-env-ref sym)))
     
    98100                                          (set-car! svec v))
    99101                                      i))
    100                   (else  i))))) i lst)) 0))
     102                  (else  i)))))) i lst)) 0))
    101103
    102104 
     
    131133               (fold
    132134                (lambda (x i)
    133                   (let-values (((df sym rhs)  (match x ((sym rhs) (values #f sym rhs))
    134                                                      ((df sym rhs) (values df sym rhs)))))
    135                     (if df (let ((svec  (solve-env-ref sym))
    136                                  (hy1  (vector-ref hy1vect i)))
    137                              (if (not (zero? hy1))
    138                                  (let* ((oldv (svec-value svec))
    139                                         (v    (+ oldv hy1)))
    140                                    (svec-value-set! svec v)
    141                                    (vector-set! z i (- hy1 (- v  oldv)))))
    142                              (fx+ 1 i))
    143                         i))) i lst)) 0)
     135                  (let ((xrec (if (> (length x) 2) x (cons #f x))))
     136                    (let ((df   (first xrec))
     137                          (sym  (second xrec))
     138                          (rhs  (third xrec)))
     139                      (if df (let ((svec  (solve-env-ref sym))
     140                                   (hy1  (vector-ref hy1vect i)))
     141                               (if (not (zero? hy1))
     142                                   (let* ((oldv (svec-value svec))
     143                                          (v    (+ oldv hy1)))
     144                                     (svec-value-set! svec v)
     145                                     (vector-set! z i (- hy1 (- v  oldv)))))
     146                               (fx+ 1 i))
     147                          i)))) i lst)) 0)
    144148
    145149            (let* ((t1  (+ t tstep tz))
     
    152156                        (print "euler: ")
    153157                        (for-each (lambda (x) (print x)) (ode-runtime 'debug))))
    154             (values it t tstep (list))))))
     158            (list it t tstep (list))))))
    155159
    156160  (if debug (begin
  • ode/trunk/examples/ccompile.scm

    r7225 r7961  
    11#!/bin/sh
    22#|  -*- Hen -*-
    3 exec sh -c "csi -qb \"$0\" -- $1 $2 $3 $4 $5 "
     3exec sh -c "csi -s \"$0\" $1 $2 $3 $4 $5 "
    44|#
    55
  • ode/trunk/examples/run.scm

    r7225 r7961  
    11#!/bin/sh
    22#|  -*- Hen -*-
    3 exec sh -c "csi -qb \"$0\" -- $1 $2 $3 $4 $5 "
     3exec sh -c "csi -s \"$0\" $1 $2 $3 $4 $5 "
    44|#
    55
  • ode/trunk/rkf45.scm

    r7924 r7961  
    152152       (for-each
    153153        (lambda (x)
    154           (let-values (((df sym rhs)  (match x ((sym rhs) (values #f sym rhs))
    155                                              ((df sym rhs) (values df sym rhs)))))
     154          (let ((xrec (if (> (length x) 2) x (cons #f x))))
     155            (let ((df   (first xrec))
     156                  (sym  (second xrec))
     157                  (rhs  (third xrec)))
    156158            (let ((v     (eval-rhs rhs))
    157159                  (svec  (solve-env-ref sym)))
     
    163165                                        (if (vector? svec)
    164166                                            (svec-value-set! svec v)
    165                                             (set-car! svec v))))))))
     167                                            (set-car! svec v)))))))))
    166168        lst))))
    167169 
     
    420422         
    421423          (begin
    422             (values it t tstep (list))))))
     424            (list it t tstep (list))))))
    423425
    424426  (if debug (print "rkf45: " (ode-runtime 'debug)))
Note: See TracChangeset for help on using the changeset viewer.