# source:git/kernel/p_Mult_q.cc@2f959e8

spielwiese
Last change on this file since 2f959e8 was 2f959e8, checked in by Oliver Wienand <wienand@…>, 16 years ago
p_Mult_q: select by domain property kutil.cc: no ntl needed anymore git-svn-id: file:///usr/local/Singular/svn/trunk@10539 2c84dea3-7e68-4137-9b89-c4e89433aadc
• Property mode set to `100644`
File size: 6.2 KB
Line
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4/***************************************************************
5 *  File:    p_Mult_q.cc
6 *  Purpose: multiplication of polynomials
7 *  Author:  obachman (Olaf Bachmann)
8 *  Created: 8/00
9 *  Version: \$Id: p_Mult_q.cc,v 1.6 2008-01-30 18:56:37 wienand Exp \$
10 *******************************************************************/
11#include "mod2.h"
12
13/***************************************************************
14 *
15 * Returns:  p * q,
16 * Destroys: if !copy then p, q
17 * Assumes: pLength(p) >= 2 pLength(q) >=2
18 ***************************************************************/
19#include "p_polys.h"
20#include "p_Procs.h"
21#include "p_Numbers.h"
22#include "kbuckets.h"
23
24#include "p_Mult_q.h"
25
26BOOLEAN pqLength(poly p, poly q, int &lp, int &lq, const int min)
27{
28  int l = 0;
29
30  do
31  {
32    if (p == NULL)
33    {
34      lp = l;
35      if (l < min)
36      {
37        if (q != NULL)
38          lq = l+1;
39        else
40          lq = l;
41        return FALSE;
42      }
43      lq = l + pLength(q);
44      return TRUE;
45    }
46    pIter(p);
47    if (q == NULL)
48    {
49      lq = l;
50      if (l < min)
51      {
52        lp = l+1;
53        return FALSE;
54      }
55      lp = l + 1 + pLength(p);
56      return TRUE;
57    }
58    pIter(q);
59    l++;
60  }
61  while (1);
62}
63
64
65static poly _p_Mult_q_Bucket(poly p, const int lp,
66                             poly q, const int lq,
67                             const int copy, const ring r)
68{
69  assume(p != NULL && pNext(p) != NULL && q != NULL && pNext(q) != NULL);
70  pAssume1(! pHaveCommonMonoms(p, q));
71  assume(lp >= 1 && lq >= 1);
72  p_Test(p, r);
73  p_Test(q, r);
74
75  poly res = pp_Mult_mm(p,q,r);     // holds initially q1*p
76#ifdef HAVE_RINGS
77  if (rField_is_Ring(currRing))
78    WarnS("_p_Mult_q_Bucket not ring testet !!!");
79#endif
80  poly qq = pNext(q);               // we iter of this
81  poly qn = pp_Mult_mm(qq, p,r);    // holds p1*qi
82  poly pp = pNext(p);               // used for Lm(qq)*pp
83  poly rr = res;                    // last monom which is surely not NULL
84  poly rn = pNext(res);             // pNext(rr)
85  number n, n1;
86
87  kBucket_pt bucket = kBucketCreate(r);
88
89  // initialize bucket
90  kBucketInit(bucket, pNext(rn), lp - 2);
91  pNext(rn) = NULL;
92
93  // now the main loop
94  Top:
95  if (rn == NULL) goto Smaller;
96  p_LmCmpAction(rn, qn, r, goto Equal, goto Greater, goto Smaller);
97
98  Greater:
99  // rn > qn, so iter
100  rr = rn;
101  pNext(rn) = kBucketExtractLm(bucket);
102  pIter(rn);
103  goto Top;
104
105  // rn < qn, append qn to rr, and compute next Lm(qq)*pp
106  Smaller:
107  pNext(rr) = qn;
108  rr = qn;
109  pIter(qn);
110  Work: // compute res + Lm(qq)*pp
111  if (rn == NULL)
112  {
113    pNext(rr) = pp_Mult_mm(pp, qq, r);
114    kBucketInit(bucket, pNext(pNext(rr)), lp - 2);
115    pNext(pNext(rr)) = NULL;
116  }
117  else
118  {
119    kBucketSetLm(bucket, rn);
120    kBucket_Plus_mm_Mult_pp(bucket, qq, pp, lp - 1);
121    pNext(rr) = kBucketExtractLm(bucket);
122  }
123
124  pIter(qq);
125  if (qq == NULL) goto Finish;
126  rn = pNext(rr);
127  goto Top;
128
129  Equal:
130  n1 = pGetCoeff(rn);
131  n = n_Add(n1, pGetCoeff(qn), r);
132  n_Delete(&n1, r);
133  if (n_IsZero(n, r))
134  {
135    n_Delete(&n, r);
136    p_LmFree(rn, r);
137  }
138  else
139  {
140    pSetCoeff0(rn, n);
141    rr = rn;
142  }
143  rn = kBucketExtractLm(bucket);
144  n_Delete(&pGetCoeff(qn),r);
145  qn = p_LmFreeAndNext(qn, r);
146  goto Work;
147
148  Finish:
149  assume(rr != NULL && pNext(rr) != NULL);
150  pNext(pNext(rr)) = kBucketClear(bucket);
151  kBucketDestroy(&bucket);
152
153  if (!copy)
154  {
155    p_Delete(&p, r);
156    p_Delete(&q, r);
157  }
158  p_Test(res, r);
159  return res;
160}
161
162#ifdef HAVE_RINGS
163static poly _p_Mult_q_Normal_ZeroDiv(poly p, poly q, const int copy, const ring r)
164{
165  assume(p != NULL && pNext(p) != NULL && q != NULL && pNext(q) != NULL);
166  pAssume1(! pHaveCommonMonoms(p, q));
167  p_Test(p, r);
168  p_Test(q, r);
169
170  poly res = pp_Mult_mm(p,q,r);     // holds initially q1*p
171  poly qq = pNext(q);               // we iter of this
172
173  while (qq != NULL)
174  {
175    res = p_Plus_mm_Mult_qq(res, qq, p, r);
176    pIter(qq);
177  }
178
179  if (!copy)
180  {
181    p_Delete(&p, r);
182    p_Delete(&q, r);
183  }
184
185  p_Test(res, r);
186
187  return res;
188}
189#endif
190
191static poly _p_Mult_q_Normal(poly p, poly q, const int copy, const ring r)
192{
193  assume(p != NULL && pNext(p) != NULL && q != NULL && pNext(q) != NULL);
194#ifdef HAVE_RINGS
195  assume(rField_is_Domain(r));
196#endif
197  pAssume1(! pHaveCommonMonoms(p, q));
198  p_Test(p, r);
199  p_Test(q, r);
200
201  poly res = pp_Mult_mm(p,q,r);     // holds initially q1*p
202  poly qq = pNext(q);               // we iter of this
203  poly qn = pp_Mult_mm(qq, p,r);    // holds p1*qi
204  poly pp = pNext(p);               // used for Lm(qq)*pp
205  poly rr = res;                    // last monom which is surely not NULL
206  poly rn = pNext(res);             // pNext(rr)
207  number n, n1;
208
209  // now the main loop
210  Top:
211  if (rn == NULL) goto Smaller;
212  p_LmCmpAction(rn, qn, r, goto Equal, goto Greater, goto Smaller);
213
214  Greater:
215  // rn > qn, so iter
216  rr = rn;
217  pIter(rn);
218  goto Top;
219
220  // rn < qn, append qn to rr, and compute next Lm(qq)*pp
221  Smaller:
222  pNext(rr) = qn;
223  rr = qn;
224  pIter(qn);
225
226  Work: // compute res + Lm(qq)*pp
227  if (rn == NULL)
228    pNext(rr) = pp_Mult_mm(pp, qq, r);
229  else
230  {
231    pNext(rr) = p_Plus_mm_Mult_qq(rn, qq, pp, r);
232  }
233
234  pIter(qq);
235  if (qq == NULL) goto Finish;
236  rn = pNext(rr);
237  goto Top;
238
239  Equal:
240  n1 = pGetCoeff(rn);
241  n = n_Add(n1, pGetCoeff(qn), r);
242  n_Delete(&n1, r);
243  if (n_IsZero(n, r))
244  {
245    n_Delete(&n, r);
246    rn = p_LmFreeAndNext(rn, r);
247  }
248  else
249  {
250    pSetCoeff0(rn, n);
251    rr = rn;
252    pIter(rn);
253  }
254  n_Delete(&pGetCoeff(qn),r);
255  qn = p_LmFreeAndNext(qn, r);
256  goto Work;
257
258  Finish:
259  if (!copy)
260  {
261    p_Delete(&p, r);
262    p_Delete(&q, r);
263  }
264  p_Test(res, r);
265  return res;
266}
267
268
269poly _p_Mult_q(poly p, poly q, const int copy, const ring r)
270{
271  int lp, lq, l;
272  poly pt;
273
274  pqLength(p, q, lp, lq, MIN_LENGTH_BUCKET);
275
276  if (lp < lq)
277  {
278    pt = p;
279    p =  q;
280    q = pt;
281    l = lp;
282    lp = lq;
283    lq = l;
284  }
285  if (lq < MIN_LENGTH_BUCKET || TEST_OPT_NOT_BUCKETS)
286#ifdef HAVE_RINGS
287    if (!rField_is_Domain(currRing))
288      return _p_Mult_q_Normal_ZeroDiv(p, q, copy, r);
289    else
290#endif
291    return _p_Mult_q_Normal(p, q, copy, r);
292  else
293  {
294    assume(lp == pLength(p));
295    assume(lq == pLength(q));
296    return _p_Mult_q_Bucket(p, lp, q, lq, copy, r);
297  }
298}
Note: See TracBrowser for help on using the repository browser.