/[emacs]/emacs/lisp/calc/calc-poly.el
ViewVC logotype

Diff of /emacs/lisp/calc/calc-poly.el

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.8 by jpb, Wed Nov 17 19:23:01 2004 UTC revision 1.9 by jpb, Thu Nov 25 05:52:38 2004 UTC
# Line 516  Line 516 
516    
517  ;;; Given an expression find all variables that are polynomial bases.  ;;; Given an expression find all variables that are polynomial bases.
518  ;;; Return list in the form '( (var1 degree1) (var2 degree2) ... ).  ;;; Return list in the form '( (var1 degree1) (var2 degree2) ... ).
519  ;;; Note dynamic scope of mpb-total-base.  
520    ;; The variable math-poly-base-total-base is local to
521    ;; math-total-polynomial-base, but is used by math-polynomial-p1,
522    ;; which is called by math-total-polynomial-base.
523    (defvar math-poly-base-total-base)
524    
525  (defun math-total-polynomial-base (expr)  (defun math-total-polynomial-base (expr)
526    (let ((mpb-total-base nil))    (let ((math-poly-base-total-base nil))
527      (math-polynomial-base expr 'math-polynomial-p1)      (math-polynomial-base expr 'math-polynomial-p1)
528      (math-sort-poly-base-list mpb-total-base)))      (math-sort-poly-base-list math-poly-base-total-base)))
529    
530    ;; The variable math-poly-base-top-expr is local to math-polynomial-base
531    ;; in calc-alg.el, but is used by math-polynomial-p1 which is called
532    ;; by math-polynomial-base.
533    (defvar math-poly-base-top-expr)
534    
535  (defun math-polynomial-p1 (subexpr)  (defun math-polynomial-p1 (subexpr)
536    (or (assoc subexpr mpb-total-base)    (or (assoc subexpr math-poly-base-total-base)
537        (memq (car subexpr) '(+ - * / neg))        (memq (car subexpr) '(+ - * / neg))
538        (and (eq (car subexpr) '^) (natnump (nth 2 subexpr)))        (and (eq (car subexpr) '^) (natnump (nth 2 subexpr)))
539        (let* ((math-poly-base-variable subexpr)        (let* ((math-poly-base-variable subexpr)
540               (exponent (math-polynomial-p mpb-top-expr subexpr)))               (exponent (math-polynomial-p math-poly-base-top-expr subexpr)))
541          (if exponent          (if exponent
542              (setq mpb-total-base (cons (list subexpr exponent)              (setq math-poly-base-total-base (cons (list subexpr exponent)
543                                         mpb-total-base)))))                                         math-poly-base-total-base)))))
544    nil)    nil)
545    
546    ;; The variable math-factored-vars is local to calcFunc-factors and
547    ;; calcFunc-factor, but is used by math-factor-expr and
548    ;; math-factor-expr-part, which are called (directly and indirectly) by
549    ;; calcFunc-factor and calcFunc-factors.
550    (defvar math-factored-vars)
551    
552    ;; The variable math-fact-expr is local to calcFunc-factors,
553    ;; calcFunc-factor and math-factor-expr, but is used by math-factor-expr-try
554    ;; and math-factor-expr-part, which are called (directly and indirectly) by
555    ;; calcFunc-factor, calcFunc-factors and math-factor-expr.
556    (defvar math-fact-expr)
557    
558    ;; The variable math-to-list is local to calcFunc-factors and
559    ;; calcFunc-factor, but is used by math-accum-factors, which is
560    ;; called (indirectly) by calcFunc-factors and calcFunc-factor.
561    (defvar math-to-list)
562    
563    (defun calcFunc-factors (math-fact-expr &optional var)
   
 (defun calcFunc-factors (expr &optional var)  
564    (let ((math-factored-vars (if var t nil))    (let ((math-factored-vars (if var t nil))
565          (math-to-list t)          (math-to-list t)
566          (calc-prefer-frac t))          (calc-prefer-frac t))
567      (or var      (or var
568          (setq var (math-polynomial-base expr)))          (setq var (math-polynomial-base math-fact-expr)))
569      (let ((res (math-factor-finish      (let ((res (math-factor-finish
570                  (or (catch 'factor (math-factor-expr-try var))                  (or (catch 'factor (math-factor-expr-try var))
571                      expr))))                      math-fact-expr))))
572        (math-simplify (if (math-vectorp res)        (math-simplify (if (math-vectorp res)
573                           res                           res
574                         (list 'vec (list 'vec res 1)))))))                         (list 'vec (list 'vec res 1)))))))
575    
576  (defun calcFunc-factor (expr &optional var)  (defun calcFunc-factor (math-fact-expr &optional var)
577    (let ((math-factored-vars nil)    (let ((math-factored-vars nil)
578          (math-to-list nil)          (math-to-list nil)
579          (calc-prefer-frac t))          (calc-prefer-frac t))
580      (math-simplify (math-factor-finish      (math-simplify (math-factor-finish
581                      (if var                      (if var
582                          (let ((math-factored-vars t))                          (let ((math-factored-vars t))
583                            (or (catch 'factor (math-factor-expr-try var)) expr))                            (or (catch 'factor (math-factor-expr-try var)) math-fact-expr))
584                        (math-factor-expr expr))))))                        (math-factor-expr math-fact-expr))))))
585    
586  (defun math-factor-finish (x)  (defun math-factor-finish (x)
587    (if (Math-primp x)    (if (Math-primp x)
# Line 571  Line 595 
595        (list 'calcFunc-Fac-Prot x)        (list 'calcFunc-Fac-Prot x)
596      x))      x))
597    
598  (defun math-factor-expr (expr)  (defun math-factor-expr (math-fact-expr)
599    (cond ((eq math-factored-vars t) expr)    (cond ((eq math-factored-vars t) math-fact-expr)
600          ((or (memq (car-safe expr) '(* / ^ neg))          ((or (memq (car-safe math-fact-expr) '(* / ^ neg))
601               (assq (car-safe expr) calc-tweak-eqn-table))               (assq (car-safe math-fact-expr) calc-tweak-eqn-table))
602           (cons (car expr) (mapcar 'math-factor-expr (cdr expr))))           (cons (car math-fact-expr) (mapcar 'math-factor-expr (cdr math-fact-expr))))
603          ((memq (car-safe expr) '(+ -))          ((memq (car-safe math-fact-expr) '(+ -))
604           (let* ((math-factored-vars math-factored-vars)           (let* ((math-factored-vars math-factored-vars)
605                  (y (catch 'factor (math-factor-expr-part expr))))                  (y (catch 'factor (math-factor-expr-part math-fact-expr))))
606             (if y             (if y
607                 (math-factor-expr y)                 (math-factor-expr y)
608               expr)))               math-fact-expr)))
609          (t expr)))          (t math-fact-expr)))
610    
611  (defun math-factor-expr-part (x)    ; uses "expr"  (defun math-factor-expr-part (x)    ; uses "expr"
612    (if (memq (car-safe x) '(+ - * / ^ neg))    (if (memq (car-safe x) '(+ - * / ^ neg))
# Line 590  Line 614 
614          (math-factor-expr-part (car x)))          (math-factor-expr-part (car x)))
615      (and (not (Math-objvecp x))      (and (not (Math-objvecp x))
616           (not (assoc x math-factored-vars))           (not (assoc x math-factored-vars))
617           (> (math-factor-contains expr x) 1)           (> (math-factor-contains math-fact-expr x) 1)
618           (setq math-factored-vars (cons (list x) math-factored-vars))           (setq math-factored-vars (cons (list x) math-factored-vars))
619           (math-factor-expr-try x))))           (math-factor-expr-try x))))
620    
621  (defun math-factor-expr-try (x)  ;; The variable math-fet-x is local to math-factor-expr-try, but is
622    (if (eq (car-safe expr) '*)  ;; used by math-factor-poly-coefs, which is called by math-factor-expr-try.
623        (let ((res1 (catch 'factor (let ((expr (nth 1 expr)))  (defvar math-fet-x)
624                                     (math-factor-expr-try x))))  
625              (res2 (catch 'factor (let ((expr (nth 2 expr)))  (defun math-factor-expr-try (math-fet-x)
626                                     (math-factor-expr-try x)))))    (if (eq (car-safe math-fact-expr) '*)
627          (let ((res1 (catch 'factor (let ((math-fact-expr (nth 1 math-fact-expr)))
628                                       (math-factor-expr-try math-fet-x))))
629                (res2 (catch 'factor (let ((math-fact-expr (nth 2 math-fact-expr)))
630                                       (math-factor-expr-try math-fet-x)))))
631          (and (or res1 res2)          (and (or res1 res2)
632               (throw 'factor (math-accum-factors (or res1 (nth 1 expr)) 1               (throw 'factor (math-accum-factors (or res1 (nth 1 math-fact-expr)) 1
633                                                  (or res2 (nth 2 expr))))))                                                  (or res2 (nth 2 math-fact-expr))))))
634      (let* ((p (math-is-polynomial expr x 30 'gen))      (let* ((p (math-is-polynomial math-fact-expr math-fet-x 30 'gen))
635             (math-poly-modulus (math-poly-modulus expr))             (math-poly-modulus (math-poly-modulus math-fact-expr))
636             res)             res)
637        (and (cdr p)        (and (cdr p)
638             (setq res (math-factor-poly-coefs p))             (setq res (math-factor-poly-coefs p))
# Line 642  Line 670 
670      (math-mul (math-pow fac pow) facs)))      (math-mul (math-pow fac pow) facs)))
671    
672  (defun math-factor-poly-coefs (p &optional square-free)    ; uses "x"  (defun math-factor-poly-coefs (p &optional square-free)    ; uses "x"
673    (let (t1 t2)    (let (t1 t2 temp)
674      (cond ((not (cdr p))      (cond ((not (cdr p))
675             (or (car p) 0))             (or (car p) 0))
676    
677            ;; Strip off multiples of x.            ;; Strip off multiples of math-fet-x.
678            ((Math-zerop (car p))            ((Math-zerop (car p))
679             (let ((z 0))             (let ((z 0))
680               (while (and p (Math-zerop (car p)))               (while (and p (Math-zerop (car p)))
# Line 654  Line 682 
682               (if (cdr p)               (if (cdr p)
683                   (setq p (math-factor-poly-coefs p square-free))                   (setq p (math-factor-poly-coefs p square-free))
684                 (setq p (math-sort-terms (math-factor-expr (car p)))))                 (setq p (math-sort-terms (math-factor-expr (car p)))))
685               (math-accum-factors x z (math-factor-protect p))))               (math-accum-factors math-fet-x z (math-factor-protect p))))
686    
687            ;; Factor out content.            ;; Factor out content.
688            ((and (not square-free)            ((and (not square-free)
# Line 665  Line 693 
693             (math-accum-factors t1 1 (math-factor-poly-coefs             (math-accum-factors t1 1 (math-factor-poly-coefs
694                                       (math-poly-div-list p t1) 'cont)))                                       (math-poly-div-list p t1) 'cont)))
695    
696            ;; Check if linear in x.            ;; Check if linear in math-fet-x.
697            ((not (cdr (cdr p)))            ((not (cdr (cdr p)))
698             (math-add (math-factor-protect             (math-add (math-factor-protect
699                        (math-sort-terms                        (math-sort-terms
700                         (math-factor-expr (car p))))                         (math-factor-expr (car p))))
701                       (math-mul x (math-factor-protect                       (math-mul math-fet-x (math-factor-protect
702                                    (math-sort-terms                                    (math-sort-terms
703                                     (math-factor-expr (nth 1 p)))))))                                     (math-factor-expr (nth 1 p)))))))
704    
# Line 683  Line 711 
711                 (setq pp (cdr pp)))                 (setq pp (cdr pp)))
712               pp)               pp)
713             (let ((res (math-rewrite             (let ((res (math-rewrite
714                         (list 'calcFunc-thecoefs x (cons 'vec p))                         (list 'calcFunc-thecoefs math-fet-x (cons 'vec p))
715                         '(var FactorRules var-FactorRules))))                         '(var FactorRules var-FactorRules))))
716               (or (and (eq (car-safe res) 'calcFunc-thefactors)               (or (and (eq (car-safe res) 'calcFunc-thefactors)
717                        (= (length res) 3)                        (= (length res) 3)
# Line 693  Line 721 
721                          (while (setq vec (cdr vec))                          (while (setq vec (cdr vec))
722                            (setq facs (math-accum-factors (car vec) 1 facs)))                            (setq facs (math-accum-factors (car vec) 1 facs)))
723                          facs))                          facs))
724                   (math-build-polynomial-expr p x))))                   (math-build-polynomial-expr p math-fet-x))))
725    
726            ;; Check if rational coefficients (i.e., not modulo a prime).            ;; Check if rational coefficients (i.e., not modulo a prime).
727            ((eq math-poly-modulus 1)            ((eq math-poly-modulus 1)
# Line 724  Line 752 
752                                             (setq scale (math-div scale den))                                             (setq scale (math-div scale den))
753                                             (math-add                                             (math-add
754                                              (math-add                                              (math-add
755                                               (math-mul den (math-pow x 2))                                               (math-mul den (math-pow math-fet-x 2))
756                                               (math-mul (math-mul coef1 den) x))                                               (math-mul (math-mul coef1 den)
757                                                           math-fet-x))
758                                              (math-mul coef0 den)))                                              (math-mul coef0 den)))
759                                         (let ((den (math-lcm-denoms coef0)))                                         (let ((den (math-lcm-denoms coef0)))
760                                           (setq scale (math-div scale den))                                           (setq scale (math-div scale den))
761                                           (math-add (math-mul den x)                                           (math-add (math-mul den math-fet-x)
762                                                     (math-mul coef0 den))))                                                     (math-mul coef0 den))))
763                                       1 expr)                                       1 expr)
764                                 roots (cdr roots))))                                 roots (cdr roots))))
# Line 738  Line 767 
767                                   (math-mul csign                                   (math-mul csign
768                                             (math-build-polynomial-expr                                             (math-build-polynomial-expr
769                                              (math-mul-list (nth 1 t1) scale)                                              (math-mul-list (nth 1 t1) scale)
770                                              x)))))                                              math-fet-x)))))
771                   (math-build-polynomial-expr p x))   ; can't factor it.                   (math-build-polynomial-expr p math-fet-x))   ; can't factor it.
772    
773               ;; Separate out the squared terms (Knuth exercise 4.6.2-34).               ;; Separate out the squared terms (Knuth exercise 4.6.2-34).
774               ;; This step also divides out the content of the polynomial.               ;; This step also divides out the content of the polynomial.

Legend:
Removed from v.1.8  
changed lines
  Added in v.1.9

savannah-hackers-public@gnu.org
ViewVC Help
Powered by ViewVC 1.1.26