home *** CD-ROM | disk | FTP | other *** search
/ Enter 2004 January / enter-2004-01.iso / files / maxima-5.9.0.exe / {app} / share / maxima / 5.9.0 / src / psolve.lisp < prev    next >
Encoding:
Text File  |  2003-02-09  |  7.1 KB  |  283 lines
  1. ;;; -*-  Mode: Lisp; Package: Maxima; Syntax: Common-Lisp; Base: 10 -*- ;;;;
  2. ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
  3. ;;;     The data in this file contains enhancments.                    ;;;;;
  4. ;;;                                                                    ;;;;;
  5. ;;;  Copyright (c) 1984,1987 by William Schelter,University of Texas   ;;;;;
  6. ;;;     All rights reserved                                            ;;;;;
  7. ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
  8. ;; ** (c) Copyright 1980 Massachusetts Institute of Technology **
  9.  
  10. (in-package "MAXIMA")
  11. (macsyma-module psolve)
  12.  
  13. (DECLARE-TOP(GENPREFIX PSO)
  14.      (SPECIAL MULT *ROOTS *FAILURES $SOLVEFACTORS))
  15.  
  16. (DECLARE-TOP(SPLITFILE SCUBIC))
  17.  
  18. (DEFMVAR FLAG4 NIL)
  19.  
  20. (DEFMFUN SOLVECUBIC (X) 
  21.        (PROG (S1 A0 A1 A2 DISCR LCOEF ADIV3 OMEGA^2 PDIV3 QDIV-2
  22.           OMEGA Y1 U y2) 
  23.          (SETQ X (CDR X))
  24.          (SETQ LCOEF (CADR X))
  25.          (SETQ ADIV3
  26.            (LIST '(MTIMES)
  27.              '((RAT) -1. 3.)
  28.              (RDIS (SETQ A2 (RATREDUCE (PTERM X 2.)
  29.                            LCOEF)))))
  30.          (SETQ A1 (RATREDUCE (PTERM X 1.) LCOEF))
  31.          (SETQ A0 (RATREDUCE (PTERM X 0.) LCOEF))
  32.          (SETQ S1 '((MTIMES)
  33.             ((RAT) 1. 2.)
  34.             $%I
  35.             ((MEXPT) 3. ((RAT) 1. 2.))))
  36.          (SETQ OMEGA (LIST '(MPLUS)
  37.                    '((RAT) -1. 2.)
  38.                    S1) 
  39.            OMEGA^2 (LIST '(MPLUS)
  40.                  '((RAT) -1. 2.)
  41.                  (LIST '(MTIMES) -1. S1)))
  42.          (SETQ PDIV3
  43.            (MEVAL* (RDIS (RATPLUS (RATTIMES A1 '(1. . 3.) T)
  44.                       (RATTIMES (RATEXPT A2 2.)
  45.                             '(-1. . 9.)
  46.                             T)))))
  47.          (AND (NOT (EQUAL PDIV3 0.)) (GO HARDER))
  48.          (SETQ 
  49.           Y1
  50.           (SIMPTIMES
  51.            (LIST
  52.         '(MTIMES)
  53.         '((RAT) 1. 3.)
  54.         (LIST '(MPLUS)
  55.                (SIMPNRT (RDIS (setq y2 (RATPLUS (RATEXPT A2 3.)
  56.                                (RATTIMES '(-27. . 1.)
  57.                                  A0
  58.                                  T))))
  59.                 3)
  60.               (LIST '(MTIMES) -1. (RDIS A2))))
  61.            1.
  62.            NIL))
  63.          (AND FLAG4 (RETURN (SOLVE3 Y1 MULT)))
  64.          (setq y2 (simpnrt (rdis (rattimes  y2 '(1. . 27.) t))
  65.                    3))
  66.          (RETURN (MAPC #'(LAMBDA (J) (SOLVE3 J MULT))
  67.                (LIST Y1
  68.                  (LIST '(MPLUS)
  69.                        (LIST '(MTIMES)
  70.                          OMEGA
  71.                          Y2)
  72.                        ADIV3)
  73.                  (LIST '(MPLUS)
  74.                        (LIST '(MTIMES)
  75.                          OMEGA^2
  76.                          Y2)
  77.                        ADIV3))))
  78.     HARDER
  79.          (SETQ 
  80.           QDIV-2
  81.           (RDIS (RATPLUS (RATTIMES (RATPLUS (RATTIMES A1 A2 T)
  82.                         (RATTIMES '(-3. . 1.)
  83.                               A0
  84.                               T))
  85.                        '(1. . 6.)
  86.                        T)
  87.                  (RATTIMES (RATEXPT A2 3.)
  88.                        '(-1. . 27.)
  89.                        T))))
  90.          (COND ((EQUAL QDIV-2 0.)
  91.             (SETQ U (SIMPNRT PDIV3 2))
  92.             (SETQ Y1 ADIV3))
  93.            (T (SETQ DISCR (SIMPLUS (LIST '(MPLUS)
  94.                          (LIST '(MEXPT)
  95.                                PDIV3
  96.                                3.)
  97.                          (LIST '(MEXPT)
  98.                                QDIV-2
  99.                                2.))
  100.                        1.
  101.                        NIL))
  102.               (COND ((EQUAL DISCR 0.)
  103.                  (SETQ U (SIMPNRT QDIV-2 3)))
  104.                 (T (SETQ DISCR (SIMPNRT DISCR 2))
  105.                    (AND (COMPLICATED DISCR)
  106.                     (SETQ DISCR (adispline DISCR)))
  107.                    (SETQ U (SIMPEXPT (LIST '(MEXPT)
  108.                          (LIST '(MPLUS)
  109.                            QDIV-2
  110.                            DISCR)
  111.                          '((RAT) 1 3)) 1 NIL))
  112.                    (AND (COMPLICATED U)
  113.                     (SETQ U (adispline U)))))))
  114.          (IF (EQUAL U 0) (MERROR "Arithmetic overflow - SOLVECUBIC"))
  115.          (OR Y1
  116.          (SETQ Y1 (SIMPLUS (LIST '(MPLUS)
  117.                      ADIV3
  118.                      U
  119.                      (LIST '(MTIMES)
  120.                            -1.
  121.                            PDIV3
  122.                            (LIST '(MEXPT)
  123.                              U
  124.                              -1.)))
  125.                    1.
  126.                    NIL)))
  127.          (RETURN
  128.           (COND (FLAG4 (SOLVE3 Y1 MULT))
  129.             (T (MAPC 
  130.             #'(LAMBDA (J) (SOLVE3 J MULT))
  131.             (LIST Y1
  132.                   (LIST '(MPLUS)
  133.                     ADIV3
  134.                     (LIST '(MTIMES) OMEGA U)
  135.                     (LIST '(MTIMES)
  136.                       -1.
  137.                       PDIV3
  138.                       OMEGA^2
  139.                       (LIST '(MEXPT)
  140.                         U
  141.                         -1.)))
  142.                   (LIST '(MPLUS)
  143.                     ADIV3
  144.                     (LIST '(MTIMES) OMEGA^2 U)
  145.                     (LIST '(MTIMES)
  146.                       -1.
  147.                       PDIV3
  148.                       OMEGA
  149.                       (LIST '(MEXPT)
  150.                         U
  151.                         -1.))))))))))
  152.  
  153. (declare-top (SPLITFILE SQUART))
  154.  
  155. (DEFMFUN SOLVEQUARTIC (X) 
  156.        (PROG (A0 A1 A2 B1 B2 B3 B0 LCOEF Z1 R TR1 TR2 D D1 E SQB3) 
  157.          (SETQ X (CDR X) LCOEF (CADR X))
  158.          (SETQ B3 (RATREDUCE (PTERM X 3.) LCOEF))
  159.          (SETQ B2 (RATREDUCE (PTERM X 2.) LCOEF))
  160.          (SETQ B1 (RATREDUCE (PTERM X 1.) LCOEF))
  161.          (SETQ B0 (RATREDUCE (PTERM X 0.) LCOEF))
  162.          (SETQ A2 (RATMINUS B2))
  163.          (SETQ A1 (RATDIF (RATTIMES B1 B3 T)
  164.                   (SETQ A0 (RATTIMES B0
  165.                          '(4. . 1.)
  166.                          T))))
  167.          (SETQ A0
  168.            (RATDIF (RATDIF (RATTIMES B2 A0 T)
  169.                    (RATTIMES (SETQ SQB3
  170.                            (RATEXPT B3 2.))
  171.                          B0
  172.                          T))
  173.                (RATEXPT B1 2.)))
  174.          (SETQ 
  175.           TR2
  176.           (SIMPLIFY (RDIS
  177.            (RATTIMES
  178.         '(1. . 4.)
  179.         (RATDIF (RATDIF (RATTIMES B3
  180.                       (RATTIMES B2
  181.                             '(4. . 1.)
  182.                             T)
  183.                       T)
  184.                 (RATTIMES '(8. . 1.) B1 T))
  185.             (RATTIMES SQB3 B3 NIL))
  186.         T))))
  187.          (SETQ Z1 (RESOLVENT A2 A1 A0))
  188.          (SETQ R
  189.            (SIMPLUS (LIST '(MPLUS)
  190.                   Z1
  191.                   (RDIS (RATDIF (RATTIMES SQB3
  192.                               '(1. . 4.)
  193.                               T)
  194.                         B2)))
  195.                 1.
  196.                 NIL))
  197.          (AND (EQUAL R 0.) (GO L0))
  198.          (SETQ R (SIMPNRT R 2))
  199.          (AND (COMPLICATED R) (SETQ R (adispline R)))
  200.          (AND (COMPLICATED TR2) (SETQ TR2 (adispline TR2)))
  201.          (SETQ TR1
  202.            (SIMPLUS (LIST '(MPLUS)
  203.                   (RDIS (RATDIF (RATTIMES SQB3
  204.                               '(1. . 2.)
  205.                               T)
  206.                         B2))
  207.                   (LIST '(MTIMES) -1. Z1))
  208.                 1.
  209.                 NIL))
  210.          (AND (COMPLICATED TR1) (SETQ TR1 (adispline TR1)))
  211.          (SETQ TR2 (DIV* TR2 R))
  212.          (GO LB1)
  213.     L0   (SETQ D1
  214.            (SIMPNRT (SIMPLIFY (LIST '(MPLUS)
  215.                   (LIST '(MEXPT) Z1 2.)
  216.                   (LIST '(MTIMES)
  217.                     -4.
  218.                     (RDIS B0))))
  219.                 2))
  220.          (SETQ TR2 (SIMPLIFY (LIST '(MTIMES) 2. D1)))
  221.          (AND (COMPLICATED TR2) (SETQ TR2 (adispline TR2)))
  222.          (SETQ TR1
  223.            (SIMPLIFY (RDIS (RATDIF (RATTIMES SQB3 '(3. . 4.) T)
  224.                  (RATTIMES B2 '(2. . 1.) T)))))
  225.          (AND (COMPLICATED TR1) (SETQ TR1 (adispline TR1)))
  226.     LB1  (SETQ D (SIMPNRT (SIMPLIFY (LIST '(MPLUS) TR1 TR2)) 2))
  227.          (SETQ E
  228.            (SIMPNRT (SIMPLIFY (LIST '(MPLUS)
  229.                   TR1
  230.                   (LIST '(MTIMES) -1. TR2)))
  231.                 2))
  232.          (SETQ D (DIV* D 2.))
  233.          (AND (COMPLICATED D) (SETQ D (adispline D)))
  234.          (SETQ E (DIV* E 2.))
  235.          (AND (COMPLICATED E) (SETQ E (adispline E)))
  236.          (SETQ A2 (RDIS (RATTIMES B3 '(-1. . 4.) T)))
  237.          (SETQ A1 (DIV* R 2.))
  238.          (SETQ Z1
  239.            (LIST (LIST '(MPLUS) A2 A1 D)
  240.              (LIST '(MPLUS)
  241.                    A2
  242.                    A1
  243.                    (LIST '(MTIMES) -1. D))
  244.              (LIST '(MPLUS)
  245.                    A2
  246.                    (LIST '(MTIMES) -1. A1)
  247.                    E)
  248.              (LIST '(MPLUS)
  249.                    A2
  250.                    (LIST '(MTIMES) -1. A1)
  251.                    (LIST '(MTIMES) -1. E))))
  252.          (RETURN (MAPC #'(LAMBDA (J) (SOLVE3 J MULT))
  253.                Z1))))
  254.  
  255. ;;; SOLVES RESOLVENT CUBIC EQUATION
  256. ;;; GENERATED FROM QUARTIC
  257.  
  258. (DEFUN RESOLVENT (A2 A1 A0) 
  259.        (PROG (*ROOTS FLAG4 *FAILURES $solvefactors)    ;undoes binding in
  260.          (SETQ FLAG4 T $solvefactors t)        ;algsys
  261.          (SOLVE (SIMPLUS (LIST '(MPLUS)
  262.                    (LIST '(MEXPT)
  263.                      'YY
  264.                      3.)
  265.                    (LIST '(MTIMES)
  266.                      (RDIS A2)
  267.                      (LIST '(MEXPT)
  268.                            'YY
  269.                            2.))
  270.                    (LIST '(MTIMES)
  271.                      (RDIS A1)
  272.                      'YY)
  273.                    (RDIS A0))
  274.                  1.
  275.                  NIL)
  276.             'YY
  277.             1.)
  278.          (COND ((zl-MEMBER 0. *ROOTS) (RETURN 0.)))
  279.          (RETURN (CADDAR (CDR (REVERSE *ROOTS))))))
  280.  
  281. #-NIL
  282. (DECLARE-TOP(UNSPECIAL MULT))
  283.