Tangent.

.


Tangent



The specificity of Lisp is, maybe, best visible on problems that require processing of data naturally represented as S-expressions. An example of such problem is "Probability that random propositional formula is tautology." In this post, the variation of classical example - symbolic differentiation - is presented: newLISP (and it is very similar in other Lisps) program that computes the tangent of the graph of the function f of single variable, defined as composition of elementary functions (+, -, sin, cos ...) in given point x0. The tangent of the function f in x0 is defined as linear function

y(x) = a·(x - x0) + b,

where a = f '(x0) and b = f(x0). The program consists of:

  • The function tangent that accepts two arguments: the function f in form of lambda list and the floating point number x0. The function tangent analyzes the code of the function f, calculates the values a and b according to mathematical definition of tangent (see above) and constructs the code of the function

    (lambda(x)(+. (*. a (-. x x0)) b))

    where symbols a, b and x0 are replaced with their respective values; that function is then returned as result of the function tangent. The implementation is very short, but it delegates main part of the work to the function d.
  • The function d accepts two arguments: formula and variable. It returns the expression obtained by symbolic differentiation of the formula for given variable. It is main function in this program, typical recursive "code walking" through formula. Simple "domain specific language" is used for descriptions of the rules  for symbolic differentiation that should be applied, for example

    (+. (*. df g) (*. f dg))

    for multiplication. The application of  rules for symbolic differentiation frequently give result that can be simplified. For instance, derivation of (*. 2 x) is

    (+. (*. 0 x) (*. 2 1))

    so function simplify is called when appropriate.
  • The function simplify accepts only one argument: symbolic expression formula. It analyzes and simplify formula, again, using typical recursive "code walking". For instance, S-expression (- x x) is simplified to 0. There are infinitely many possible rules for simplification; the function performs only few of these. The function uses eval at one place: if some expression contains only operators and constants, but not variables then the easisest way to simplify it is to eval it. (One might speculate that simplify is generalized eval.) Simplification is, actually, not necessary for computation of the tangent. However, its use in the context of symbolic differentiation is reasonable.
Few operators from my library are used; function expand// i.e. parallel expand, fexpr println= for convenient printing and floating point arithmetic operators +., -., *. and /., synonymous for built-in add, sub, mul and div. These functions are not essential.

The graph of the complicated function and tangent on the curve suggests that program, generally, works.

(setf f0 (lambda(x)
           (+. (sin (*. 12 x))
               (cos (*. 32 x))
               (tan (*. x 1.4))
               (asin x)
               (acos  x)
               (atan x)
               (*. x (cos (/. 7 x)))
               (sqrt (/. 9 x))
               (pow x x)
               (*. x (sinh x))
               (*. x (cosh x))
               (asinh x)
               (sin (acosh (+. x 1)))
               (atanh x))))
(setf x0 0.4)
(println= (tangent f0 x0))
; (tangent f0 x0)=(lambda (x) (+. (*. -21.491 (-. x 0.4)) 10.252))


graph


Finally, the code:
(setf [print.supressed] true [println.supressed] true)
(load "http://www.instprog.com/Instprog.default-library.lsp")
(setf [print.supressed] nil [println.supressed] nil)

(define (tangent f x0)
   (letn((variable           (first (first f)))
         (expression         (last f))
         (derived-expression (d expression variable))
         
         (a                  (eval (expand derived-expression
                                           (list (list variable x0)))))
         (b                  (f x0)))
        
       (expand// '(lambda(x)(+. (*. a (-. x x0)) b))
                  'a 'b 'x0)))

(define (d formula variable)
  (simplify 
    (cond 
      ((= formula variable) 1)
      ((atom? formula) 0)
      ((list? formula)
       (letn((operator (first formula))
             (operands (rest formula))
             (lexpand 
               (lambda(expr)
                 (letn((flatexpr (flat expr))
                       (f  (if (find 'f  flatexpr)(operands 0)))
                       (df (if (find 'df flatexpr)(d f variable)))
                       (g  (if (find 'g  flatexpr)(operands 1)))
                       (dg (if (find 'dg flatexpr)(d g variable))))
                   (expand// expr 'f 'df 'g 'dg)))))

         (case operator 

           (+. (cons '+. (map (lambda(op)(d op variable)) operands)))
           (-. (cons '-. (map (lambda(op)(d op variable)) operands)))

           (*. (case (length operands)
                 (1    (lexpand 'df))
                 (2    (lexpand '(+. (*. df g) (*. f dg))))
                 (true (d (list '*. (first operands)
                                    (cons '*. (rest operands)))
                           variable))))

           (/. (case (length operands)
                 (1    (d (list '/. 1 (first operands)) variable)) 
                 (2    (lexpand '(/. (-. (*. df g) (*. f dg)) (*. g g))))
                 (true (d (list '/.  (first operands)
                                     (cons '*. (rest operands)))
                       variable))))

           (pow (d (lexpand '(exp (*. g (log f)))) variable))
           (exp (lexpand '(*. f df)))
           (log (if (= (length operands) 1)
                    (lexpand '(/. df f))
                    (d (lexpand '(/. (log f) (log g))) variable)))

           (sqrt  (lexpand '(*. 0.5 df (/. 1 (sqrt f)))))
           (sin   (lexpand '(*. (cos f) df)))
           (cos   (lexpand '(*. (-. (sin f)) df)))
           (tan   (lexpand '(/. df (pow (cos f) 2))))
            (asin (lexpand '(/. df (sqrt (-. 1 (*. f f))))))
           (acos  (lexpand '(-. (/. df (sqrt (-. 1 (*. f f)))))))
           (atan  (lexpand '(/. df (+. 1 (*. f f)))))
           (sinh  (lexpand '(*. (cosh f) df)))
           (cosh  (lexpand '(*. (sinh f) df)))
           (tanh  (lexpand '(*. (-. 1 (pow (tanh f) 2)) df)))
           (asinh (lexpand '(/. df (sqrt (+.(*. f f) 1)))))
           (acosh (lexpand '(/. df (sqrt (-. (*. f f) 1)))))
           (atanh (lexpand '(/. df (-. 1 (*. f f)))))

           ))))))

(define (simplify formula)
  (cond 
    ((atom? formula) formula)
    ((list? formula)
     (letn((operator (first formula))
           (operands (map simplify (rest formula)))
           (formula (cons operator operands)))

       (cond 
       
         ; if all operands are constants, then 
         ; simplified formula is evaluated formula
   
         ((for-all number? operands)(eval formula))
   
         ; (*. x), (+. x) => x
         
         ((and (or (= operator '*.) (= operator '+.)) 
               (= (length operands) 1))
           (first operands))
           
         ; (*. ... 0 ...) => 0  
           
         ((and (= operator '*.) (find 0 operands)) 0)
         
         ; (*. ... 1 ...) => (*. ...)
         
         ((and (= operator '*.) (find 1 operands))
          (simplify (clean (curry = 1) formula)))
         
         ; (+. ... 0 ...) => 0
         
         ((and (= operator '+.) (find 0 operands))
          (simplify (clean zero? formula)))
          
         ; (-. (-. ...)) => ...
         
         ((match '(-. (-. ?)) formula)
           (last (last formula)))
           
         ; (-. minuend ...)
         
         ((and (= operator '-.) (> (length operands) 1))
               
          (letn((minuend (first operands))
                (subtrahends (rest operands))
                (subtrahend  (simplify (cons '+. subtrahends))))
                
            (cond ((zero? minuend)        (simplify (list '-. 
                                                          subtrahend)))
                  ((zero? subtrahend)     minuend)
                  ((= minuend subtrahend) 0)
                  (true                   (list '-. minuend 
                                                    subtrahend)))))
                     
         ; (/. (/. ...))
               
         ((match '(/. (/. ?)) formula) (last (last formula)))
         
         ; (/. dividend ...)
         
         ((and (= operator '/.) (> (length operands) 1))
          (letn((dividend (first operands))
                (divisors (rest operands))
                (divisor  (simplify (cons '*. divisors))))
            
            (cond ((zero? dividend)     0)
                  ((= divisor 1)        dividend)
                  ((= divisor -1)       (simplify (list '-. dividend)))
                  ((= dividend divisor) 1)
                  (true                 (list '/. dividend divisor)))))
          
         (true formula))))))
No comments:

Conflation of Subtraction and Additive Inverse in Lisp.



Conflation of Subtraction and Additive Inverse in Lisp

In mathematics, the symbol - is used as name of two different operations: subtraction and additive inverse. For instance, that symbol has different meaning in the expressions (-3)+7 and (6-3)+7. There is no ambiguity. The edge case of subtraction

c - b1 - b2 - ... - bn

for n = 0, i.e. subtraction with zero subtrahends, is c, which is clearly different than additive inverse of c: - c. However, in prefix notation use of the same symbol for two different operators will cause ambiguity so (- c) can be interpreted as both additive inverse of c and the edge case of subtraction with zero subtrahends. Designers of almost all Lisp dialects have chosen to keep the same symbol and implement both operations. The decision is based on the number of operands: if there is only one operand, operator - behaves as additive inverse; if there is two or more operands - behaves as subtraction.

This decision rarely causes any problems with "hand written" expressions because humans automatically, without much thinking, reduce edge case of subtraction to minuend, so, we'll not see the edge case of subtraction in any mathematical formulas. However, if expressions are processed by programs the reduction of the edge case doesn't happen automatically. For instance, the program that simplifies arithmetic expression could contain the function that deletes zeroes from subtrahends. That function should accept expression as

(- c b1 0 b2 0 0 0)

as argument and return expression

(- c b1 b2).

However, simple implementation like this one in Newlisp 
(define (simplify E)
  (let ((operator (first E))
        (minuend (first (rest E)))
        (subtrahends (rest (rest E))))
       (setf new-subtrahends (clean zero? subtrahends))
       (append (list operator)
               (list minuend)
               new-subtrahends)))
will not work correctly, because, for instance, (- 3 0) will be reduced into (- 3) and these two expressions have not the same value. Explicit treatment of the edge case is required:
(define (simplify E)
  (let ((operator (first E))
        (minuend (first (rest E)))
        (subtrahends (rest (rest E))))
       (setf new-subtrahends (clean zero? subtrahends))
       (if (empty? new-subtrahends)
           minuend
           (append (list operator)
                   (list minuend)
                   new-subtrahends))))
The second definition is significantly (~ 15%) larger and less consistent: one can get rid of subtrahends equal to zero only if he removes whole subtraction. On the other hand, edge cases for some other operations like addition are well supported in all Lisp dialects: not only that (+ a 0) can be safely simplified to (+ a), but even (+ 0 0) can be simplified to (+).

Hence, merging of the two operators into one could be design mistake. It would be, arguably, better to
  • define new operator for additive inverse and
  • extend definition of subtraction on special case with zero subtrahends.
Of course, S-expressions would look even less similar to arithmetic expressions, however, clarity and consistency over convenience is, I think, the element of the philosophy of Lisp.

It is even more surprising that this merging is generalized: division operator, / which is not used as divisional inverse in mathematical expressions is used on the way analogous to - in almost all Lisp dialects.

It is hard to expect change in existing Lisp dialects because lot of existing code uses "conflated operator", although I think that such changes should be, gradually, done: the past is finite and future is infinite. However, I think that for future Lisp dialects it would be better to keep the difference between subtraction (division) and additive (multiplicative) inverse.
2 comments:
Subscribe to: Posts (Atom)