source: git/kernel/clapconv.cc @ 1629ab

spielwiese
Last change on this file since 1629ab was 1629ab, checked in by Hans Schönemann <hannes@…>, 15 years ago
*hannes: cleanup longalg git-svn-id: file:///usr/local/Singular/svn/trunk@12071 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 12.9 KB
Line 
1// emacs edit mode for this file is -*- C++ -*-
2/****************************************
3*  Computer Algebra System SINGULAR     *
4****************************************/
5// $Id: clapconv.cc,v 1.17 2009-08-14 17:14:29 Singular Exp $
6/*
7* ABSTRACT: convert data between Singular and factory
8*/
9
10
11#include "mod2.h"
12#ifdef HAVE_FACTORY
13#include "omalloc.h"
14#include "structs.h"
15#define SI_DONT_HAVE_GLOBAL_VARS
16#include "clapconv.h"
17#include "numbers.h"
18#include "modulop.h"
19#include "longalg.h"
20#include "polys.h"
21#include "febase.h"
22#include "ring.h"
23#include "sbuckets.h"
24#include "ffields.h"
25void out_cf(char *s1,const CanonicalForm &f,char *s2);
26
27static void convRec( const CanonicalForm & f, int * exp, poly & result );
28
29static void convRecAlg( const CanonicalForm & f, int * exp, napoly & result );
30
31static void conv_RecPP ( const CanonicalForm & f, int * exp, sBucket_pt result, ring r );
32
33static void convRecTrP ( const CanonicalForm & f, int * exp, poly & result, int offs, const ring r );
34
35static void convRecGFGF ( const CanonicalForm & f, int * exp, poly & result );
36
37static number convFactoryNSingAN( const CanonicalForm &f);
38
39CanonicalForm convSingNFactoryN( number n )
40{
41  CanonicalForm term;
42  /* does only work for Zp, Q */
43  if ( getCharacteristic() != 0 )
44  {
45    term = npInt( n );
46  }
47  else
48  {
49    if ( SR_HDL(n) & SR_INT )
50    {
51      term = SR_TO_INT(n);
52    }
53    else
54    {
55      if ( n->s == 3 )
56      {
57        MP_INT dummy;
58        mpz_init_set( &dummy, &(n->z) );
59        term = make_cf( dummy );
60      }
61      else
62      {
63        // assume s==0 or s==1
64        MP_INT num, den;
65        On(SW_RATIONAL);
66        mpz_init_set( &num, &(n->z) );
67        mpz_init_set( &den, &(n->n) );
68        term = make_cf( num, den, ( n->s != 1 ));
69      }
70    }
71  }
72  return term;
73}
74
75number convFactoryNSingN( const CanonicalForm & n)
76{
77  if (n.isImm())
78    return nInit(n.intval());
79  else
80  {
81    number z=(number)omAllocBin(rnumber_bin);
82#if defined(LDEBUG)
83    z->debug=123456;
84#endif
85    z->z = gmp_numerator( n );
86    if ( n.den().isOne() )
87      z->s = 3;
88    else
89    {
90      z->n = gmp_denominator( n );
91      z->s = 0;
92    }
93    return z;
94  }
95}
96
97poly convFactoryPSingP ( const CanonicalForm & f, const ring r )
98{
99//    cerr << " f = " << f << endl;
100  int n = rVar(r)+1;
101  /* ASSERT( level( f ) <= pVariables, "illegal number of variables" ); */
102  int * exp = (int*)omAlloc0(n*sizeof(int));
103  sBucket_pt result_bucket=sBucketCreate(r);
104  conv_RecPP( f, exp, result_bucket, r );
105  poly result; int dummy;
106  sBucketDestroyMerge(result_bucket,&result,&dummy);
107  omFreeSize((ADDRESS)exp,n*sizeof(int));
108  return result;
109}
110
111static void conv_RecPP ( const CanonicalForm & f, int * exp, sBucket_pt result, ring r )
112{
113  if (f.isZero())
114    return;
115  if ( ! f.inCoeffDomain() )
116  {
117    int l = f.level();
118    for ( CFIterator i = f; i.hasTerms(); i++ )
119    {
120      exp[l] = i.exp();
121      conv_RecPP( i.coeff(), exp, result, r );
122    }
123    exp[l] = 0;
124  }
125  else
126  {
127    poly term = p_Init(r);
128    pNext( term ) = NULL;
129    for ( int i = 1; i <= r->N; i++ )
130      p_SetExp( term, i, exp[i], r);
131    //if (rRing_has_Comp(r)) p_SetComp(term, 0, r); // done by p_Init
132    if ( f.isImm() )
133      pGetCoeff( term ) = n_Init( f.intval(), r );
134    else
135    {
136      number z=(number)omAllocBin(rnumber_bin);
137#if defined(LDEBUG)
138      z->debug=123456;
139#endif
140      z->z = gmp_numerator( f );
141      if ( f.den().isOne() )
142        z->s = 3;
143      else
144      {
145        z->n = gmp_denominator( f );
146        z->s = 0;
147        nlNormalize(z);
148      }
149      pGetCoeff( term ) = z;
150    }
151    p_Setm( term, r );
152    if ( n_IsZero(pGetCoeff(term), r) )
153    {
154      p_Delete(&term,r);
155    }
156    else
157    {
158      sBucket_Merge_p(result,term,1);
159    }
160  }
161}
162
163
164CanonicalForm convSingPFactoryP( poly p, const ring r )
165{
166  CanonicalForm result = 0;
167  int e, n = rVar(r);
168
169  while ( p!=NULL)
170  {
171    CanonicalForm term;
172    /* does only work for finite fields */
173    if ( rField_is_Zp(r) )
174    {
175      term = npInt( pGetCoeff( p ) );
176    }
177    else if (rField_is_Q(r))
178    {
179      if ( SR_HDL(pGetCoeff( p )) & SR_INT )
180        term = SR_TO_INT( pGetCoeff( p ) );
181      else
182      {
183        if ( pGetCoeff( p )->s == 3 )
184        {
185          MP_INT dummy;
186          mpz_init_set( &dummy, &(pGetCoeff( p )->z) );
187          term = make_cf( dummy );
188        }
189        else
190        {
191          // assume s==0 or s==1
192          MP_INT num, den;
193          On(SW_RATIONAL);
194          mpz_init_set( &num, &(pGetCoeff( p )->z) );
195          mpz_init_set( &den, &(pGetCoeff( p )->n) );
196          term = make_cf( num, den, ( pGetCoeff( p )->s != 1 ));
197        }
198      }
199    }
200    else
201    {
202      WerrorS("conversion error");
203      return result;
204    }
205    for ( int i = n; i >0; i-- )
206    {
207      if ( (e = p_GetExp( p, i, r)) != 0 )
208        term *= power( Variable( i ), e );
209    }
210    result += term;
211    pIter( p );
212  }
213  return result;
214}
215
216CanonicalForm convSingAPFactoryAP ( poly p , const Variable & a)
217{
218  CanonicalForm result = 0;
219  int e, n = pVariables;
220  int off=rPar(currRing);
221
222  On(SW_RATIONAL);
223  while ( p!=NULL)
224  {
225    CanonicalForm term=convSingAFactoryA(((lnumber)pGetCoeff(p))->z,a);
226    for ( int i = 1; i <= n; i++ )
227    {
228      if ( (e = pGetExp( p, i )) != 0 )
229        term *= power( Variable( i + off), e );
230    }
231    result += term;
232    pIter( p );
233  }
234  return result;
235}
236
237static void
238convRecAP_R ( const CanonicalForm & f, int * exp, poly & result, int par_start, int var_start) ;
239
240poly convFactoryAPSingAP_R ( const CanonicalForm & f, int par_start, int var_start )
241{
242  int n = pVariables+1+rPar(currRing);
243  int * exp = (int *)omAlloc0(n*sizeof(int));
244  poly result = NULL;
245  convRecAP_R( f, exp, result,par_start, var_start );
246  omFreeSize((ADDRESS)exp,n*sizeof(int));
247  return result;
248}
249
250poly convFactoryAPSingAP ( const CanonicalForm & f )
251{
252  return convFactoryAPSingAP_R(f,0,rPar(currRing));
253}
254
255static void convRecAP_R ( const CanonicalForm & f, int * exp, poly & result, int par_start, int var_start )
256{
257  if (f.isZero())
258    return;
259  if ( ! f.inCoeffDomain() )
260  {
261    int l = f.level();
262    for ( CFIterator i = f; i.hasTerms(); i++ )
263    {
264      exp[l] = i.exp();
265      convRecAP_R( i.coeff(), exp, result, par_start, var_start );
266    }
267    exp[l] = 0;
268  }
269  else
270  {
271    napoly z=(napoly)convFactoryASingA( f );
272    if (z!=NULL)
273    {
274      poly term = pInit();
275      pNext( term ) = NULL;
276      int i;
277      for ( i = 1; i <= pVariables; i++ )
278        pSetExp( term, i , exp[i+var_start]);
279      //if (rRing_has_Comp(currRing->algring)) p_SetComp(term, 0, currRing->algring); // done by pInit
280      if (par_start==0)
281      {
282        for ( i = 1; i <= var_start; i++ )
283        //z->e[i-1]+=exp[i];
284          napAddExp(z,i,exp[i]);
285      }
286      else
287      {
288        for ( i = par_start+1; i <= var_start+rPar(currRing); i++ )
289        //z->e[i-1]+=exp[i];
290          napAddExp(z,i,exp[i-par_start]);
291      }
292      pGetCoeff(term)=(number)omAlloc0Bin(rnumber_bin);
293      ((lnumber)pGetCoeff(term))->z=z;
294      pSetm( term );
295      result = pAdd( result, term );
296    }
297  }
298}
299
300CanonicalForm convSingAFactoryA ( napoly p , const Variable & a )
301{
302  CanonicalForm result = 0;
303  int e;
304
305  while ( p!=NULL )
306  {
307    CanonicalForm term;
308    /* does only work for finite fields:Z/p */
309    if ( rField_is_Zp_a() )
310    {
311      term = npInt( napGetCoeff( p ) );
312    }
313    else
314    {
315      if ( SR_HDL(napGetCoeff( p )) & SR_INT )
316        term = SR_TO_INT(napGetCoeff( p )) ;
317      else
318      {
319        if ( napGetCoeff( p )->s == 3 )
320        {
321          MP_INT dummy;
322          mpz_init_set( &dummy, &(napGetCoeff( p )->z) );
323          term = make_cf( dummy );
324        }
325        else
326        {
327          // assume s==0 or s==1
328          MP_INT num, den;
329          On(SW_RATIONAL);
330          mpz_init_set( &num, &(napGetCoeff( p )->z) );
331          mpz_init_set( &den, &(napGetCoeff( p )->n) );
332          term = make_cf( num, den, ( napGetCoeff( p )->s != 1 ));
333        }
334      }
335    }
336    if ( (e = napGetExp( p, 1 )) != 0 )
337      term *= power( a , e );
338    result += term;
339    p = pNext( p );
340  }
341  return result;
342}
343
344static number convFactoryNSingAN( const CanonicalForm &f)
345{
346  if ( f.isImm() )
347    return n_Init( f.intval(), currRing->algring );
348  else
349  {
350    number z=(number)omAllocBin(rnumber_bin);
351#if defined(LDEBUG)
352    z->debug=123456;
353#endif
354    z->z = gmp_numerator( f );
355    if ( f.den().isOne() )
356    {
357      z->s = 3;
358    }
359    else
360    {
361      z->n = gmp_denominator( f );
362      z->s = 0;
363      nlNormalize(z);
364    }
365    #ifdef LDEBUG
366    nlDBTest(z,__FILE__,__LINE__);
367    #endif
368    return z;
369  }
370}
371
372napoly convFactoryASingA ( const CanonicalForm & f )
373{
374  napoly a=NULL;
375  napoly t;
376  for( CFIterator i=f; i.hasTerms(); i++)
377  {
378    t=napNew();
379    // pNext( t ) = NULL; //already done by napNew
380    napGetCoeff(t)=convFactoryNSingAN( i.coeff() );
381    if (n_IsZero(napGetCoeff(t),currRing->algring))
382    {
383      p_Delete(&t,nacRing);
384    }
385    else
386    {
387      napSetExp(t,1,i.exp());
388      a=napAdd(a,t);
389    }
390  }
391  if (a!=NULL)
392  {
393    if (naMinimalPoly!=NULL)
394    {
395      if (napGetExp(a,1) >= napGetExp(naMinimalPoly,1))
396        a = napRemainder( a, naMinimalPoly);
397    }
398  }
399  return a;
400}
401
402CanonicalForm convSingTrPFactoryP ( poly p, const ring r )
403{
404  CanonicalForm result = 0;
405  int e, n = rVar(r);
406  int offs = rPar(r);
407
408  while ( p!=NULL )
409  {
410    n_Normalize(pGetCoeff(p),r);
411    CanonicalForm term=convSingPFactoryP(((lnumber)pGetCoeff(p))->z,r->algring);
412
413    if ((((lnumber)pGetCoeff(p))->n!=NULL)
414    && (!errorreported))
415    {
416      WerrorS("conversion error: denominator!= 1");
417    }
418
419    for ( int i = n; i > 0; i-- )
420    {
421      if ( (e = p_GetExp( p, i,r )) != 0 )
422        term = term * power( Variable( i + offs ), e );
423    }
424    result += term;
425    p = pNext( p );
426  }
427  return result;
428}
429
430poly convFactoryPSingTrP ( const CanonicalForm & f, const ring r )
431{
432  int n = rVar(r)+1;
433  int * exp = (int*)omAlloc0(n*sizeof(int));
434  poly result = NULL;
435  convRecTrP( f, exp, result , rPar(r), r );
436  omFreeSize((ADDRESS)exp,n*sizeof(int));
437  return result;
438}
439
440static void
441convRecTrP ( const CanonicalForm & f, int * exp, poly & result , int offs, const ring r)
442{
443  if (f.isZero())
444    return;
445  if ( f.level() > offs )
446  {
447    int l = f.level();
448    for ( CFIterator i = f; i.hasTerms(); i++ )
449    {
450      exp[l-offs] = i.exp();
451      convRecTrP( i.coeff(), exp, result, offs, r );
452    }
453    exp[l-offs] = 0;
454  }
455  else
456  {
457    poly term = p_Init(r);
458    pNext( term ) = NULL;
459    for ( int i = rVar(r); i>0; i-- )
460      p_SetExp( term, i ,exp[i], r);
461    //if (rRing_has_Comp(currRing)) p_SetComp(term, 0, currRing); // done by pInit
462    pGetCoeff(term)=(number)omAlloc0Bin(rnumber_bin);
463    ((lnumber)pGetCoeff(term))->z=convFactoryPSingP( f, r->algring );
464    p_Setm( term,r );
465    result = p_Add_q( result, term,r );
466  }
467}
468
469number   nlChineseRemainder(number *x, number *q,int rl)
470// elemenst in the array are x[0..(rl-1)], q[0..(rl-1)]
471{
472#ifdef HAVE_FACTORY
473  CFArray X(rl), Q(rl);
474  int i;
475  for(i=rl-1;i>=0;i--)
476  {
477    X[i]=CanonicalForm(nlInt(x[i]));
478    Q[i]=CanonicalForm(nlInt(q[i]));
479  }
480  CanonicalForm xnew,qnew;
481  chineseRemainder(X,Q,xnew,qnew);
482  number n=convFactoryNSingN(xnew);
483  number p=convFactoryNSingN(qnew);
484  number p2=nlIntDiv(p,nlInit(2));
485  if (nlGreater(n,p2))
486  {
487     number n2=nlSub(n,p);
488     nlDelete(&n,currRing);
489     n=n2;
490  }
491  nlDelete(&p,currRing);
492  nlDelete(&p2,currRing);
493  return n;
494#else
495  WerrorS("not implemented");
496  return nlInit(0);
497#endif
498}
499
500CanonicalForm
501convSingGFFactoryGF( poly p )
502{
503  CanonicalForm result=CanonicalForm(0);
504  int e, n = pVariables;
505
506  while ( p != NULL )
507  {
508    CanonicalForm term;
509    term = make_cf_from_gf( (int)(long)pGetCoeff( p ) );
510    //int * A=(int *)&term;
511    //Print("term=%x, == 0 ?: %d\n",*A, term.isZero());
512    for ( int i = 1; i <= n; i++ )
513    {
514      if ( (e = pGetExp( p, i )) != 0 )
515         term *= power( Variable( i ), e );
516    }
517    result += term;
518    p = pNext( p );
519  }
520  return result;
521}
522
523poly
524convFactoryGFSingGF ( const CanonicalForm & f )
525{
526//    cerr << " f = " << f << endl;
527  int n = pVariables+1;
528  /* ASSERT( level( f ) <= pVariables, "illegal number of variables" ); */
529  int * exp = (int*)omAlloc0(n*sizeof(int));
530  poly result = NULL;
531  convRecGFGF( f, exp, result );
532  omFreeSize((ADDRESS)exp,n*sizeof(int));
533  return result;
534}
535
536static void
537convRecGFGF ( const CanonicalForm & f, int * exp, poly & result )
538{
539  if (f.isZero())
540    return;
541  if ( ! f.inCoeffDomain() )
542  {
543    int l = f.level();
544    for ( CFIterator i = f; i.hasTerms(); i++ )
545    {
546      exp[l] = i.exp();
547      convRecGFGF( i.coeff(), exp, result );
548    }
549    exp[l] = 0;
550  }
551  else
552  {
553    poly term = pInit();
554    pNext( term ) = NULL;
555    for ( int i = 1; i <= pVariables; i++ )
556      pSetExp( term, i, exp[i]);
557    //if (rRing_has_Comp(currRing)) p_SetComp(term, 0, currRing); // done by pInit
558    pGetCoeff( term ) = (number) gf_value (f);
559    pSetm( term );
560    result = pAdd( result, term );
561  }
562}
563
564int convFactoryISingI( const CanonicalForm & f)
565{
566  if (!f.isImm()) WerrorS("int overflow in det");
567  return f.intval();
568}
569#endif
Note: See TracBrowser for help on using the repository browser.