← Back to branch summary

~ubuntu-branches/debian/squeeze/maxima/squeeze

« back to all changes in this revision

Viewing changes to share/affine/aquotient.lisp

  • Committer: Bazaar Package Importer
  • Author(s): Camm Maguire
  • Date: 2006-10-18 14:52:42 UTC
  • mto: (1.1.5 upstream)
  • mto: This revision was merged to the branch mainline in revision 4.
  • Revision ID: james.westby@ubuntu.com-20061018145242-vzyrm5hmxr8kiosf
Import upstream version 5.10.0

Show diffs side-by-side

added added

removed removed

Lines of Context:
share/affine/aquotient.lisp
5
5
;;;     All rights reserved                                            ;;;;;
6
6
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
7
7
 
8
 
(in-package "MAXIMA")
 
8
(in-package :maxima)
9
9
;leading coefficient of f
10
10
 
11
11
 
175
175
             do 
176
176
         (setq terms-g g-orig)
177
177
         (prog ()
178
 
            FIRST-PRODUCT
 
178
            first-product
179
179
               (cond ((null terms-g)(return nil)))
180
180
               (setq prod-exp (f+ f-exp (term-deg terms-g)))
181
181
               (setq prod-cof (afp-times(term-cof terms-g) f-cof))
182
 
               (cond ((pzerop prod-cof) (setq terms-g (cddr terms-g)) (GO FIRST-PRODUCT))
 
182
               (cond ((pzerop prod-cof) (setq terms-g (cddr terms-g)) (go first-product))
183
183
                     ((or (null answ) (> prod-exp (term-deg answ)))
184
184
                      (setq answ (cons prod-exp (cons prod-cof answ)))
185
 
                      (setq tail (cdr answ))                   (GO NEXT-PRODUCT))
 
185
                      (setq tail (cdr answ))                   (go next-product))
186
186
                     ((eql prod-exp (term-deg answ))
187
187
                      (setq prod-cof (pplus (term-cof answ) prod-cof))
188
188
                      (cond ((pzerop prod-cof) (setq answ (cddr answ))
189
 
                             (setq terms-g (cddr terms-g))          (GO FIRST-PRODUCT))
 
189
                             (setq terms-g (cddr terms-g))          (go first-product))
190
190
                            (t (setf (term-cof answ) prod-cof)
191
 
                               (setq tail (cdr answ))        (GO NEXT-PRODUCT)))))
 
191
                               (setq tail (cdr answ))        (go next-product)))))
192
192
               ;;below here assume answ not empty and (term-deg answ)
193
193
               ;;greater any possible future prod-exp (until next
194
194
               ;;f-exp,f-cof)  Once below this point we stay below,
196
196
               ;;coefficient whose corresponding degree is definitely
197
197
               ;;higher than any prod-exp to be encountered.
198
198
               (setq tail (cdr answ))
199
 
            TAIL-CERTAIN 
 
199
            tail-certain 
200
200
               (cond ((and (cdr tail)(> (second tail) prod-exp))
201
 
                      (setq tail (cddr tail))               (GO TAIL-CERTAIN)))
 
201
                      (setq tail (cddr tail))               (go tail-certain)))
202
202
               (cond ((or (null (cdr tail))(< (second tail) prod-exp))
203
203
                      (setf (cdr tail) (cons prod-exp (cons prod-cof (cdr tail))))
204
 
                      (setq tail (cddr tail))                (GO NEXT-PRODUCT)))
 
204
                      (setq tail (cddr tail))                (go next-product)))
205
205
               (cond ((pzerop (setq prod-cof (pplus (third tail) prod-cof)))
206
206
                      (setf (cdr tail) (cdddr tail)))
207
207
                     (t (setf (third tail) prod-cof) (setq  tail (cddr tail))))
208
 
            NEXT-PRODUCT 
 
208
            next-product 
209
209
               (setq terms-g (cddr terms-g))
210
210
               (cond ((null terms-g) (return nil)))
211
211
               (setq prod-exp (f+ f-exp (car terms-g)))
212
212
               (setq prod-cof (afp-times (second terms-g) f-cof))
213
 
               (cond ((pzerop prod-cof)                    (GO NEXT-PRODUCT)))
214
 
               (GO TAIL-CERTAIN)))))
 
213
               (cond ((pzerop prod-cof)                    (go next-product)))
 
214
               (go tail-certain)))))
215
215
  answ)
216
216
 
217
217
;;;destructive:
1197
1197
         (cond ((zerop g)(values f 1 0))
1198
1198
               (t (values 1 0  (crecip g)))))
1199
1199
        ((not (eql (p-var g) (p-var f)))(ferror 'not-function-of-one-variable))
1200
 
        (t(multiple-value-bind (quot zl-REM creqd)
 
1200
        (t(multiple-value-bind (quot zl-rem creqd)
1201
1201
              (afp-pseudo-quotient f g)
1202
1202
            creqd ;;ignore
1203
1203
            (multiple-value-bind (gcd a b)
1204
 
                (recursive-ideal-gcd1 g zl-REM)
 
1204
                (recursive-ideal-gcd1 g zl-rem)
1205
1205
              (values gcd b (pdifference a (ptimes quot b))))))))
1206
1206
 
1207
1207
 
1245
1245
                 when (equal u v)
1246
1246
                   collecting 1 into cof
1247
1247
                 else
1248
 
                 when (and (LISTP v) (MEMBER u  (cdr v) :test 'equal))
 
1248
                 when (and (listp v) (member u  (cdr v) :test 'equal))
1249
1249
                  collecting
1250
1250
                    (cond ((eql (length v) 3)
1251
1251
                           (sloop for vv in (cdr v )
1374
1374
          do (return k)
1375
1375
       do (setq u u2 v v2 k (f1+ k)))
1376
1376
  (prog ()
1377
 
     B2
 
1377
     b2
1378
1378
        (cond ((oddp u) (setq tt (- v)))
1379
1379
              (t(setq tt (ash u -1))))
1380
 
     B3B4
 
1380
     b3b4
1381
1381
        (sloop  do (setq t2 (ash tt -1))
1382
1382
               when  (eql (ash t2 1) tt)
1383
1383
               do (setq tt t2)
1399
1399
  (cond (row (cond ((< (array-total-size row)
1400
1400
                       (length poly))(adjust-array row (f+ 10 (poly-length poly))
1401
1401
                                                   :fill-pointer (fill-pointer  row)))))
1402
 
        (t (setq row (MAKE-ARRAY (f+ 10 (setq leng (poly-length poly))) :fill-pointer 0 :adjustable t))))
 
1402
        (t (setq row (make-array (f+ 10 (setq leng (poly-length poly))) :fill-pointer 0 :adjustable t))))
1403
1403
  (cond ((numberp poly)(vector-push  0 row)
1404
1404
         (vector-push  poly row))
1405
1405
        (t (sloop for u in (cdr poly)
1481
1481
              (t (return nil)))
1482
1482
     tail-certain ;;add in term to tail
1483
1483
        (cond ((and (cdr tail)(> (second tail) prod-exp))
1484
 
               (setq tail (cddr tail))               (GO TAIL-CERTAIN)))
 
1484
               (setq tail (cddr tail))               (go tail-certain)))
1485
1485
        (cond ((or (null (cdr tail))(< (second tail) prod-exp))
1486
1486
               (setf (cdr tail) (cons prod-exp (cons prod-cof (cdr tail))))
1487
 
               (setq tail (cddr tail))                (GO NEXT-DOUBLE-PRODUCT)))
 
1487
               (setq tail (cddr tail))                (go next-double-product)))
1488
1488
        (cond ((pzerop (setq prod-cof (pplus (third tail) prod-cof)))
1489
1489
               (setf (cdr tail) (cdddr tail)))
1490
1490
              (t (setf (third tail) prod-cof) (setq  tail (cddr tail))))
1509
1509
 
1510
1510
 
1511
1511
(defmacro def-test (f1 f2)
1512
 
  `(defun ,(intern (format nil "TEST-~A" f1)) (&rest rest-args)
 
1512
  `(defun ,(intern (format nil "~A-~A" '#:test f1)) (&rest rest-args)
1513
1513
     (let (empty (*print-level* 2)(*print-length* 3) ansa ansb)
1514
1514
       
1515
1515
       
1611
1611
 
1612
1612
  ;;want to find sol'ns of v^p-v=0 (mod u)
1613
1613
(defun berlekamp-set-up-and-reduce-matrix ( u p &aux rows powers estimated-size  sp)
1614
 
  (setq rows (MAKE-ARRAY (p-deg u) :fill-pointer (p-deg u)))
 
1614
  (setq rows (make-array (p-deg u) :fill-pointer (p-deg u)))
1615
1615
  
1616
1616
  (let ((modulus p)(tellratlist (list u)))
1617
1617
    (working-modulo (list u)
1629
1629
          summing (or (and (atom vv) 0) (length vv)) into count
1630
1630
          finally (setq estimated-size (f+ 10 (quotient count 5))))
1631
1631
    (sloop for i below (fill-pointer rows)
1632
 
          do (setf (aref rows i) (MAKE-ARRAY estimated-size :fill-pointer 0 :adjustable t)))
 
1632
          do (setf (aref rows i) (make-array estimated-size :fill-pointer 0 :adjustable t)))
1633
1633
    ;;putting the entries in the sparse matrix.  Each polynomial is a column.
1634
1634
    (sloop for vv  in (cdr powers)
1635
1635
          for i from 1
1695
1695
                        (setq factor-list (cons tem
1696
1696
                                                (cons (afp-quotient u-fact
1697
1697
                                                                    tem)
1698
 
                                                      (zl-DELETE u-fact factor-list))))
 
1698
                                                      (zl-delete u-fact factor-list))))
1699
1699
                        when (eql (length factor-list) number-of-factors)
1700
1700
                          do (return-from sue factor-list)))))))
1701
1701
 
1850
1850
          (sloop named sue for i from 1
1851
1851
                do
1852
1852
            (sloop  do (setq tt (generate-t-for-one-degree-factors u deg p i))
1853
 
                   unless (or  (numberp tt) (zl-MEMBER tt used-tt))
 
1853
                   unless (or  (numberp tt) (zl-member tt used-tt))
1854
1854
                     do (push tt used-tt) (return tt))
1855
1855
            (working-modulo (list u)
1856
1856
              (setq tt  (afp-expt tt pow) ))
1864
1864
                            (push (afp-make-monic tem) answ))
1865
1865
;                          ((numberp tem))
1866
1866
                           (t (push tem facts)))
1867
 
                     (setq facts(zl-DELETE v facts))))
 
1867
                     (setq facts(zl-delete v facts))))
1868
1868
                  when (null facts) do (return-from sue answ)))))
1869
1869
  (cond ((memq 1 answ) (fsignal 'bad)))
1870
1870
  (iassert (eql (p-deg  u)
1929
1929
 
1930
1930
 
1931
1931
(defun gen-afp-times (&rest lis)
1932
 
  (setq lis (zl-DELETE 1 (copy-list lis)))
 
1932
  (setq lis (zl-delete 1 (copy-list lis)))
1933
1933
  (cond ((null lis) 1)
1934
1934
        ((null (cdr lis))(car lis))
1935
1935
        (t (afp-times (car lis) (cond ((cddr lis)
1957
1957
                                    prime up-to-size))
1958
1958
            (cons (car lift) (hensel-lift-list (second lift) (cdr factor-list) prime up-to-size)))))
1959
1959
 
1960
 
(defun hensel-lift1 (product ve we prev-modulus prime  &optional a b &aux dif h kk   quot zl-REM creqd new-modulus gcd)
 
1960
(defun hensel-lift1 (product ve we prev-modulus prime  &optional a b &aux dif h kk   quot zl-rem creqd new-modulus gcd)
1961
1961
  "lifts u=ve*we mod (p^e) to u=ve+1*we+1 mod (p^e+1) with ve=ve+1 and we=we+1 mod (p^e+1)
1962
1962
  and deg(ve+1)<=deg(ve) deg(we+1)<=deg(we)  and returns the list of  ve+1 and we+1"
1963
1963
;  (declare (values (list ve+1 we+1)))
1973
1973
    (setq h (ptimes b dif))
1974
1974
    (setq kk (ptimes a dif))
1975
1975
    (let ((modulus new-modulus))
1976
 
      (multiple-value-setq ( quot zl-REM creqd)
 
1976
      (multiple-value-setq ( quot zl-rem creqd)
1977
1977
                           (afp-pseudo-quotient h ve)))
1978
 
    (setq h zl-REM)
 
1978
    (setq h zl-rem)
1979
1979
    (setq kk (pplus kk (ptimes quot we)))
1980
1980
    (list (pplus ve h) (pplus we kk))))))
1981
1981
   ;(let ((modulus (expt  prime (f1+ e)))) (values (pplus ve h) (pplus we kk))))
2056
2056
        (cond ((> (* d 2) (length factors))
2057
2057
               (push pol answ) (return answ)))
2058
2058
        (sloop for v in (sub-lists d factors)
2059
 
              when (zl-MEMBER v tried)
 
2059
              when (zl-member v tried)
2060
2060
                do nil
2061
2061
              else
2062
2062
                do (push v tried)
2066
2066
                          (setq pol quot)
2067
2067
                          (push prod answ)
2068
2068
                          (sloop for vv in v
2069
 
                                do (setq factors (zl-DELETE vv factors 1)))
 
2069
                                do (setq factors (zl-delete vv factors 1)))
2070
2070
                          (go look-for-factors)))
2071
2071
              finally (incf d)
2072
2072
                      (go look-for-factors))))
2111
2111
        (t  (setq case0 (express-x^i f g 0))
2112
2112
            (setq a (ptimes (car case0) (setq mon (list (car f) k 1))))
2113
2113
            (setq b (ptimes (second case0) mon))
2114
 
            (multiple-value-setq (quot zl-REM cre) (afp-pseudo-quotient a g))
2115
 
            (setq a zl-REM)
 
2114
            (multiple-value-setq (quot zl-rem cre) (afp-pseudo-quotient a g))
 
2115
            (setq a zl-rem)
2116
2116
            (setq b (pplus b
2117
2117
                           (ptimes quot f)))
2118
2118
            ;;(iassert (equal mon (afp-plus (afp-times f a) (afp-times g b))))
2132
2132
(defun wr-lift (u fi gi up-to-k big-mod  point &aux (modulus big-mod)(main-var (p-var u))tem fi+1 gi+1 v f0 g0 varl w
2133
2133
                aw bw )
2134
2134
  (setq v (psublis (subs-translate-sublis point) 1  u))
2135
 
  (setq varl (zl-DELETE (car u) (list-variables u)))
 
2135
  (setq varl (zl-delete (car u) (list-variables u)))
2136
2136
  (setq f0 fi )
2137
2137
  (setq g0 gi )
2138
2138
  (sloop for i from 1 to up-to-k