Changeset 15425 in project


Ignore:
Timestamp:
08/13/09 08:26:36 (10 years ago)
Author:
Ivan Raikov
Message:

added symbolic differentiation to nemo-utils

File:
1 edited

Legend:

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

    r14732 r15425  
    2424             if-convert let-enum let-elim let-lift
    2525             s+ sw+ sl\ nl spaces ppf
    26              transitions-graph state-lineqs)
     26             transitions-graph state-lineqs
     27             differentiate simplify )
    2728
    2829 (import scheme chicken data-structures srfi-1 srfi-13)
     
    3031 (require-extension matchable strictly-pretty
    3132                    varsubst digraph nemo-core)
    32 
    33 (declare (lambda-lift))
    34 
    3533
    3634
     
    214212                              lineqs)))
    215213    (list n lineqs1)))
     214
     215;;    `(+ - * / pow neg abs atan asin acos sin cos exp ln
     216;;      sqrt tan cosh sinh tanh hypot gamma lgamma log10 log2 log1p ldexp cube
     217;;      > < <= >= = and or round ceiling floor max min
     218;;      fpvector-ref))
     219
     220(define LOG10E 0.434294481903252)
     221(define LOG2E  1.44269504088896)
     222
     223(define (differentiate x t)
     224  (cond ((number? t)  0.0)
     225        ((symbol? t)  (if (equal? x t) 1.0 0.0))
     226        (else (match t
     227                (('neg u)  `(neg ,(differentiate x u)))
     228
     229                (('+ u v)  `(+ ,(differentiate x u) ,(differentiate x v)))
     230                (('- u v)  `(- ,(differentiate x u) ,(differentiate x v)))
     231
     232                (('* (and u (? number?)) v)    `(* ,u ,(differentiate x v)))
     233                (('* v (and u (? number?)))    `(* ,u ,(differentiate x v)))
     234
     235                (('* u v)     `(+ (* ,(differentiate x u) ,v)
     236                                  (* ,u ,(differentiate x v))))
     237
     238                (('/ u v)     `(/ (- (* ,(differentiate x u) ,v)
     239                                     (* ,u ,(differentiate x v)))
     240                                  (pow ,v 2.0)))
     241
     242                (('cube u)     (differentiate x `(pow ,u 3.0)))
     243
     244                (('pow u n)    (chain x u `(* ,n (pow ,u (- ,n 1.0)))))
     245               
     246                (('sqrt u)     (chain x u `(/ 1.0 (* 2.0 (sqrt ,u)))))
     247
     248                (('exp u)      (chain x u `(exp ,u)))
     249
     250                (('log a u)    (chain x u `(/ 1.0 ,u)))
     251
     252                (('log10 u)    (chain x u `(* ,LOG10E (/ ,(differentiate x u) ,u))))
     253       
     254                (('log2 u)     (chain x u `(* ,LOG2E (/ ,(differentiate x u) ,u))))
     255
     256                (('log1p u)    (differentiate x `(log (+ 1.0 ,u))))
     257
     258                (('ldexp u n)  (differentiate x `(* ,u ,(expt 2 n))))
     259       
     260                (('sin u)      (chain x u `(cos ,u)))
     261               
     262                (('cos u)      (chain x u `(neg (sin ,u))))
     263
     264                (('tan u)      (differentiate x `(* (sin ,u) (/ 1.0 (cos ,u)))))
     265               
     266                (('asin u)     (chain x u `(/ 1.0 (sqrt (- 1.0 (pow ,u 2.0))))))
     267                                             
     268                (('acos u)     (chain x u `(/ (neg 1.0) (sqrt (- 1.0 (pow ,u 2.0))))))
     269                                             
     270                (('atan u)     (chain x u `(/ 1.0 (+ 1.0 (pow ,u 2.0)))))
     271                                             
     272                (('sinh u)     (differentiate x `(/ (- (exp ,u) (exp (neg ,u))) 2.0)))
     273
     274                (('cosh u)     (differentiate x `(/ (+ (exp ,u) (exp (neg ,u))) 2.0)))
     275               
     276                (('tanh u)     (differentiate x `(/ (sinh ,u) (cosh ,u))))
     277
     278                (else          #f)))))
     279
     280(define (chain x t u)
     281  (if (symbol? t) u
     282      `(* ,(differentiate x t) ,u)))
     283       
     284(define (simplify t)
     285  (match t
     286         (('+ 0.0 t1)  t1)
     287         (('+ t1 0.0)  t1)
     288         (('+ t1 ('neg t2))  `(- ,t1 ,t2))
     289         
     290         (('- 0.0 t1)  `(neg ,t1))
     291         (('- t1 0.0)  t1)
     292         (('neg ('neg t1))  t1)
     293
     294         (('* 0.0 t1)  0.0)
     295         (('* t1 0.0)  0.0)
     296         (('* 1.0 t1)  t1)
     297         (('* t1 1.0)  t1)
     298         (('* ('neg t1) ('neg t2))  `(* ,t1 ,t2))
     299
     300         (else t)))
     301
    216302)
Note: See TracChangeset for help on using the changeset viewer.