source: git/libpolys/polys/clapconv.cc @ 75e126

spielwiese
Last change on this file since 75e126 was 6105f4f, checked in by Hans Schoenemann <hannes@…>, 9 years ago
transext.cc: NUM/DEN do not include type cast any more
  • Property mode set to 100644
File size: 10.5 KB
Line 
1// emacs edit mode for this file is -*- C++ -*-
2/****************************************
3*  Computer Algebra System SINGULAR     *
4****************************************/
5/*
6* ABSTRACT: convert data between Singular and factory
7*/
8
9
10
11
12
13#include <misc/auxiliary.h>
14#include <omalloc/omalloc.h>
15
16#include <factory/factory.h>
17
18#include <coeffs/coeffs.h>
19
20#include <coeffs/longrat.h> // snumber is necessary
21
22#include <polys/monomials/p_polys.h>
23#include <polys/sbuckets.h>
24#include <polys/clapconv.h>
25
26#include "simpleideals.h"
27
28#define TRANSEXT_PRIVATES
29#include <polys/ext_fields/transext.h>
30
31void out_cf(const char *s1,const CanonicalForm &f,const char *s2);
32
33static void conv_RecPP ( const CanonicalForm & f, int * exp, sBucket_pt result, ring r );
34
35static void convRecTrP ( const CanonicalForm & f, int * exp, poly & result, int offs, const ring r );
36
37//static void convRecGFGF ( const CanonicalForm & f, int * exp, poly & result );
38
39static number convFactoryNSingAN( const CanonicalForm &f, const ring r);
40
41poly convFactoryPSingP ( const CanonicalForm & f, const ring r )
42{
43  int n = rVar(r)+1;
44  /* ASSERT( level( f ) <= pVariables, "illegal number of variables" ); */
45  int * exp = (int*)omAlloc0(n*sizeof(int));
46  sBucket_pt result_bucket=sBucketCreate(r);
47  conv_RecPP( f, exp, result_bucket, r );
48  poly result; int dummy;
49  sBucketDestroyMerge(result_bucket,&result,&dummy);
50  omFreeSize((ADDRESS)exp,n*sizeof(int));
51  return result;
52}
53
54static void conv_RecPP ( const CanonicalForm & f, int * exp, sBucket_pt result, ring r )
55{
56  if (f.isZero())
57    return;
58  if ( ! f.inCoeffDomain() )
59  {
60    int l = f.level();
61    for ( CFIterator i = f; i.hasTerms(); i++ )
62    {
63      exp[l] = i.exp();
64      conv_RecPP( i.coeff(), exp, result, r );
65    }
66    exp[l] = 0;
67  }
68  else
69  {
70    poly term = p_Init(r);
71    pNext( term ) = NULL;
72    for ( int i = 1; i <= r->N; i++ )
73      p_SetExp( term, i, exp[i], r);
74    pGetCoeff( term )=r->cf->convFactoryNSingN(f, r->cf);
75    p_Setm( term, r );
76    if ( n_IsZero(pGetCoeff(term), r->cf) )
77    {
78      p_Delete(&term,r);
79    }
80    else
81    {
82      sBucket_Merge_p(result,term,1);
83    }
84  }
85}
86
87
88CanonicalForm convSingPFactoryP( poly p, const ring r )
89{
90  CanonicalForm result = 0;
91  int e, n = rVar(r);
92  BOOLEAN setChar=TRUE;
93
94  p=pReverse(p);
95  poly op=p;
96  while ( p!=NULL )
97  {
98    CanonicalForm term;
99    term=r->cf->convSingNFactoryN(pGetCoeff( p ),setChar, r->cf);
100    if (errorreported) break;
101    setChar=FALSE;
102    for ( int i = n; i >0; i-- )
103    {
104      if ( (e = p_GetExp( p, i, r)) != 0 )
105        term *= power( Variable( i ), e );
106    }
107    result += term;
108    pIter( p );
109 }
110 op=pReverse(op);
111 return result;
112}
113
114int convFactoryISingI( const CanonicalForm & f)
115{
116  if (!f.isImm()) WerrorS("int overflow in det");
117  return f.intval();
118}
119
120CanonicalForm convSingAPFactoryAP ( poly p , const Variable & a, const ring r)
121{
122  CanonicalForm result = 0;
123  int e, n = r-> N;
124  int off=rPar(r);
125
126  if (!rField_is_Zp_a(r))
127    On(SW_RATIONAL);
128  while ( p!=NULL)
129  {
130    CanonicalForm term=convSingAFactoryA(((poly)p_GetCoeff(p, r->cf->extRing)),a, r);
131    for ( int i = 1; i <= n; i++ )
132    {
133      if ( (e = p_GetExp( p, i, r )) != 0 )
134        term *= power( Variable( i + off), e );
135    }
136    result += term;
137    pIter( p );
138  }
139  return result;
140}
141
142static void
143convRecAP_R ( const CanonicalForm & f, int * exp, poly & result, int par_start, int var_start, const ring r) ;
144
145poly convFactoryAPSingAP_R ( const CanonicalForm & f, int par_start, int var_start, const ring r )
146{
147  int n = rVar(r)+rPar(r)+1;
148  int * exp = (int *)omAlloc0(n*sizeof(int));
149  poly result = NULL;
150  convRecAP_R( f, exp, result,par_start, var_start, r );
151  omFreeSize((ADDRESS)exp,n*sizeof(int));
152  return result;
153}
154
155poly convFactoryAPSingAP ( const CanonicalForm & f, const ring r )
156{
157  return convFactoryAPSingAP_R(f,0,rPar(r),r);
158}
159
160static void convRecAP_R ( const CanonicalForm & f, int * exp, poly & result, int par_start, int var_start, const ring r )
161{
162  if (f.isZero())
163    return;
164  if ( ! f.inCoeffDomain() )
165  {
166    int l = f.level();
167    for ( CFIterator i = f; i.hasTerms(); i++ )
168    {
169      exp[l] = i.exp();
170      convRecAP_R( i.coeff(), exp, result, par_start, var_start, r);
171    }
172    exp[l] = 0;
173  }
174  else
175  {
176    poly z=(poly)convFactoryASingA( f,r );
177    if (z!=NULL)
178    {
179      poly term = p_Init(r);
180      pNext( term ) = NULL;
181      int i;
182      for ( i = rVar(r); i>0 ; i-- )
183        p_SetExp( term, i , exp[i+var_start],r);
184      //if (rRing_has_Comp(currRing->extRing)) p_SetComp(term, 0, currRing->extRing); // done by pInit
185      if (par_start==0)
186      {
187        for ( i = 1; i <= var_start; i++ )
188        //z->e[i-1]+=exp[i];
189          p_AddExp(z,i,exp[i],r->cf->extRing);
190      }
191      else
192      {
193        for ( i = par_start+1; i <= var_start+rPar(r); i++ )
194        //z->e[i-1]+=exp[i];
195          p_AddExp(z,i,exp[i-par_start],r->cf->extRing);
196      }
197      pGetCoeff(term)=(number)ALLOC0_RNUMBER();
198      p_GetCoeff(term, r->cf->extRing)=(number) z;
199      p_Setm( term,r );
200      result = p_Add_q( result, term, r );
201    }
202  }
203}
204
205CanonicalForm convSingAFactoryA ( poly p , const Variable & a, const ring r )
206{
207  CanonicalForm result = 0;
208  int e;
209
210  while ( p!=NULL )
211  {
212    CanonicalForm term;
213    if ( rField_is_Zp_a(r) )
214    {
215      term = n_Int( p_GetCoeff( p, r->cf->extRing ), r->cf->extRing->cf );
216    }
217    else
218    {
219      if ( SR_HDL(p_GetCoeff( p, r->cf->extRing )) & SR_INT )
220        term = SR_TO_INT(p_GetCoeff( p, r->cf->extRing )) ;
221      else
222      {
223        if ( p_GetCoeff( p, r->cf->extRing )->s == 3 )
224        {
225          mpz_t dummy;
226          mpz_init_set( dummy, (p_GetCoeff( p,r->cf->extRing )->z) );
227          term = make_cf( dummy );
228        }
229        else
230        {
231          // assume s==0 or s==1
232          mpz_t num, den;
233          On(SW_RATIONAL);
234          mpz_init_set( num, (p_GetCoeff( p, r->cf->extRing )->z) );
235          mpz_init_set( den, (p_GetCoeff( p, r->cf->extRing )->n) );
236          term = make_cf( num, den, ( p_GetCoeff( p, r->cf->extRing )->s != 1 ));
237        }
238      }
239    }
240    if ( (e = p_GetExp( p, 1, r->cf->extRing )) != 0 )
241      term *= power( a , e );
242    result += term;
243    p = pNext( p );
244  }
245  return result;
246}
247
248static number convFactoryNSingAN( const CanonicalForm &f, const ring r)
249{
250  assume (r != NULL);
251  assume (r->cf != NULL);
252  assume (r->cf->extRing != NULL);
253  // r->cf->extRing->cf has to be Q or Z/p (the supported types of factory)
254  return n_convFactoryNSingN( f, r->cf->extRing->cf );
255}
256
257poly convFactoryASingA ( const CanonicalForm & f, const ring r )
258{
259  poly a=NULL;
260  poly t;
261  for( CFIterator i=f; i.hasTerms(); i++)
262  {
263    t= p_Init (r->cf->extRing);
264    p_GetCoeff(t, r->cf->extRing)= convFactoryNSingAN( i.coeff(), r );
265    if (n_IsZero(p_GetCoeff(t,r->cf->extRing),r->cf->extRing->cf))
266    {
267      p_Delete(&t,r->cf->extRing);
268    }
269    else
270    {
271      p_SetExp(t,1,i.exp(),r->cf->extRing);
272      a=p_Add_q(a,t,r->cf->extRing);
273    }
274  }
275  if (a!=NULL)
276  {
277    if( r->cf->extRing != NULL )
278      if (r->cf->extRing->qideal->m[0]!=NULL)
279      {
280        poly l=r->cf->extRing->qideal->m[0];
281        if (p_GetExp(a,1,r->cf->extRing) >= p_GetExp(l,1,r->cf->extRing))
282          a = p_PolyDiv (a, l, FALSE, r->cf->extRing); // ???
283      }
284  }
285  return a;
286}
287
288CanonicalForm convSingTrPFactoryP ( poly p, const ring r )
289{
290  CanonicalForm result = 0;
291  int e, n = rVar(r);
292  int offs = rPar(r);
293
294  while ( p!=NULL )
295  {
296    n_Normalize(p_GetCoeff(p, r), r->cf);
297
298    // test if denominator is constant
299    if (!p_IsConstantPoly(DEN ((fraction)p_GetCoeff (p,r)),r->cf->extRing) && !errorreported)
300      WerrorS("conversion error: denominator!= 1");
301
302    CanonicalForm term=convSingPFactoryP(NUM ((fraction)p_GetCoeff(p, r)),r->cf->extRing);
303
304    // if denominator is not NULL it should be a constant at this point
305    if (DEN ((fraction)p_GetCoeff(p,r)) != NULL)
306    {
307      CanonicalForm den= convSingPFactoryP(DEN ((fraction)p_GetCoeff(p, r)),r->cf->extRing);
308      if (rChar (r) == 0)
309        On (SW_RATIONAL);
310      term /= den;
311    }
312
313    for ( int i = n; i > 0; i-- )
314    {
315      if ( (e = p_GetExp( p, i,r )) != 0 )
316        term = term * power( Variable( i + offs ), e );
317    }
318    result += term;
319    p = pNext( p );
320  }
321  return result;
322}
323
324poly convFactoryPSingTrP ( const CanonicalForm & f, const ring r )
325{
326  int n = rVar(r)+1;
327  int * exp = (int*)omAlloc0(n*sizeof(int));
328  poly result = NULL;
329  convRecTrP( f, exp, result , rPar(r), r );
330  omFreeSize((ADDRESS)exp,n*sizeof(int));
331  return result;
332}
333
334static void
335convRecTrP ( const CanonicalForm & f, int * exp, poly & result , int offs, const ring r)
336{
337  if (f.isZero())
338    return;
339  if ( f.level() > offs )
340  {
341    int l = f.level();
342    for ( CFIterator i = f; i.hasTerms(); i++ )
343    {
344      exp[l-offs] = i.exp();
345      convRecTrP( i.coeff(), exp, result, offs, r );
346    }
347    exp[l-offs] = 0;
348  }
349  else
350  {
351    poly term = p_Init(r);
352    pNext( term ) = NULL;
353    for ( int i = rVar(r); i>0; i-- )
354      p_SetExp( term, i ,exp[i], r);
355    //if (rRing_has_Comp(currRing)) p_SetComp(term, 0, currRing); // done by pInit
356    pGetCoeff(term)=ntInit(convFactoryPSingP( f, r->cf->extRing ), r->cf);
357    p_Setm( term,r );
358    result = p_Add_q( result, term,r );
359  }
360}
361
362#if 0
363CanonicalForm
364convSingGFFactoryGF( poly p )
365{
366  CanonicalForm result=CanonicalForm(0);
367  int e, n = pVariables;
368
369  while ( p != NULL )
370  {
371    CanonicalForm term;
372    term = make_cf_from_gf( (int)(long)pGetCoeff( p ) );
373    //int * A=(int *)&term;
374    //Print("term=%x, == 0 ?: %d\n",*A, term.isZero());
375    for ( int i = 1; i <= n; i++ )
376    {
377      if ( (e = pGetExp( p, i )) != 0 )
378         term *= power( Variable( i ), e );
379    }
380    result += term;
381    p = pNext( p );
382  }
383  return result;
384}
385
386poly
387convFactoryGFSingGF ( const CanonicalForm & f )
388{
389//    cerr << " f = " << f << endl;
390  int n = pVariables+1;
391  /* ASSERT( level( f ) <= pVariables, "illegal number of variables" ); */
392  int * exp = (int*)omAlloc0(n*sizeof(int));
393  poly result = NULL;
394  convRecGFGF( f, exp, result );
395  omFreeSize((ADDRESS)exp,n*sizeof(int));
396  return result;
397}
398
399static void
400convRecGFGF ( const CanonicalForm & f, int * exp, poly & result )
401{
402  if (f.isZero())
403    return;
404  if ( ! f.inCoeffDomain() )
405  {
406    int l = f.level();
407    for ( CFIterator i = f; i.hasTerms(); i++ )
408    {
409      exp[l] = i.exp();
410      convRecGFGF( i.coeff(), exp, result );
411    }
412    exp[l] = 0;
413  }
414  else
415  {
416    poly term = pInit();
417    pNext( term ) = NULL;
418    for ( int i = 1; i <= pVariables; i++ )
419      pSetExp( term, i, exp[i]);
420    //if (rRing_has_Comp(currRing)) p_SetComp(term, 0, currRing); // done by pInit
421    pGetCoeff( term ) = (number) gf_value (f);
422    pSetm( term );
423    result = pAdd( result, term );
424  }
425}
426
427#endif
Note: See TracBrowser for help on using the repository browser.