source: git/kernel/clapconv.cc @ a41f3aa

spielwiese
Last change on this file since a41f3aa was 5d68a6, checked in by Hans Schönemann <hannes@…>, 15 years ago
*hannes: factory conversion:GF git-svn-id: file:///usr/local/Singular/svn/trunk@11253 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 15.6 KB
Line 
1// emacs edit mode for this file is -*- C++ -*-
2/****************************************
3*  Computer Algebra System SINGULAR     *
4****************************************/
5// $Id: clapconv.cc,v 1.12 2008-12-17 15:07:46 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 convRecPTr ( const CanonicalForm & f, int * exp, napoly & result );
34
35static void convRecTrP ( const CanonicalForm & f, int * exp, poly & result, int offs );
36
37static void convRecGFGF ( const CanonicalForm & f, int * exp, poly & result );
38
39static number convFactoryNSingAN( const CanonicalForm &f);
40
41CanonicalForm convSingNFactoryN( number n )
42{
43  CanonicalForm term;
44  /* does only work for Zp, Q */
45  if ( getCharacteristic() != 0 )
46  {
47    term = npInt( n );
48  }
49  else
50  {
51    if ( SR_HDL(n) & SR_INT )
52    {
53      term = SR_TO_INT(n);
54    }
55    else
56    {
57      if ( n->s == 3 )
58      {
59        MP_INT dummy;
60        mpz_init_set( &dummy, &(n->z) );
61        term = make_cf( dummy );
62      }
63      else
64      {
65        // assume s==0 or s==1
66        MP_INT num, den;
67        On(SW_RATIONAL);
68        mpz_init_set( &num, &(n->z) );
69        mpz_init_set( &den, &(n->n) );
70        term = make_cf( num, den, ( n->s != 1 ));
71      }
72    }
73  }
74  return term;
75}
76
77number convFactoryNSingN( const CanonicalForm & n)
78{
79  if (n.isImm())
80    return nInit(n.intval());
81  else
82  {
83    number z=(number)omAllocBin(rnumber_bin);
84#if defined(LDEBUG)
85    z->debug=123456;
86#endif
87    z->z = gmp_numerator( n );
88    if ( n.den().isOne() )
89      z->s = 3;
90    else
91    {
92      z->n = gmp_denominator( n );
93      z->s = 0;
94    }
95    return z;
96  }
97}
98
99CanonicalForm conv_SingPFactoryP( poly p, ring r )
100{
101  CanonicalForm result = 0;
102  int e, n = r->N;
103  assume( rPar(r)==0);
104
105  if ( getCharacteristic() != 0 )
106  {
107    /* does only work for finite fields Z/p*/
108    while ( p != NULL )
109    {
110      CanonicalForm term = npInt( pGetCoeff( p ) );
111      for ( int i = n; i >0; i-- )
112      {
113        if ( (e = p_GetExp( p, i, r )) != 0 )
114           term *= power( Variable( i ), e );
115      }
116      result += term;
117      pIter( p );
118    }
119  }
120  else
121  {
122    while ( p != NULL )
123    {
124      CanonicalForm term;
125      if ( SR_HDL(pGetCoeff( p )) & SR_INT )
126      {
127        term = SR_TO_INT( pGetCoeff( p ) );
128      }
129      else
130      {
131        if ( pGetCoeff( p )->s == 3 )
132        {
133          MP_INT dummy;
134          mpz_init_set( &dummy, &(pGetCoeff( p )->z) );
135          term = make_cf( dummy );
136        }
137        else
138        {
139          // assume s==1 or s==0
140          MP_INT num, den;
141          On(SW_RATIONAL);
142          mpz_init_set( &num, &(pGetCoeff( p )->z) );
143          mpz_init_set( &den, &(pGetCoeff( p )->n) );
144          term = make_cf( num, den, ( pGetCoeff( p )->s != 1 ));
145        }
146      }
147      for ( int i = n; i >0; i-- )
148      {
149        if ( (e = p_GetExp( p, i, r )) != 0 )
150           term *= power( Variable( i ), e );
151      }
152      result += term;
153      pIter( p );
154    }
155  }
156  return result;
157}
158
159poly conv_FactoryPSingP ( const CanonicalForm & f, ring r )
160{
161//    cerr << " f = " << f << endl;
162  int n = r->N+1;
163  /* ASSERT( level( f ) <= pVariables, "illegal number of variables" ); */
164  int * exp = (int*)omAlloc0(n*sizeof(int));
165  sBucket_pt result_bucket=sBucketCreate(r);
166  conv_RecPP( f, exp, result_bucket, r );
167  poly result; int dummy;
168  sBucketDestroyMerge(result_bucket,&result,&dummy);
169  omFreeSize((ADDRESS)exp,n*sizeof(int));
170  return result;
171}
172
173static void conv_RecPP ( const CanonicalForm & f, int * exp, sBucket_pt result, ring r )
174{
175  if (f.isZero())
176    return;
177  if ( ! f.inCoeffDomain() )
178  {
179    int l = f.level();
180    for ( CFIterator i = f; i.hasTerms(); i++ )
181    {
182      exp[l] = i.exp();
183      conv_RecPP( i.coeff(), exp, result, r );
184    }
185    exp[l] = 0;
186  }
187  else
188  {
189    poly term = p_Init(r);
190    pNext( term ) = NULL;
191    for ( int i = 1; i <= r->N; i++ )
192      p_SetExp( term, i, exp[i], r);
193    //if (rRing_has_Comp(r)) p_SetComp(term, 0, r); // done by p_Init
194    if ( f.isImm() )
195      pGetCoeff( term ) = n_Init( f.intval(), r );
196    else
197    {
198      number z=(number)omAllocBin(rnumber_bin);
199#if defined(LDEBUG)
200      z->debug=123456;
201#endif
202      z->z = gmp_numerator( f );
203      if ( f.den().isOne() )
204        z->s = 3;
205      else
206      {
207        z->n = gmp_denominator( f );
208        z->s = 0;
209        nlNormalize(z);
210      }
211      pGetCoeff( term ) = z;
212    }
213    p_Setm( term, r );
214    if ( n_IsZero(pGetCoeff(term), r) )
215    {
216      p_Delete(&term,r);
217    }
218    else
219    {
220      sBucket_Merge_p(result,term,1);
221    }
222  }
223}
224
225
226CanonicalForm convSingTrFactoryP( napoly p )
227{
228  CanonicalForm result = 0;
229  int e, n = napVariables;
230
231  while ( p!=NULL)
232  {
233    CanonicalForm term;
234    /* does only work for finite fields */
235    if ( getCharacteristic() != 0 )
236    {
237      term = npInt( napGetCoeff( p ) );
238    }
239    else
240    {
241      if ( SR_HDL(napGetCoeff( p )) & SR_INT )
242        term = SR_TO_INT( napGetCoeff( p ) );
243      else
244      {
245        if ( napGetCoeff( p )->s == 3 )
246        {
247          MP_INT dummy;
248          mpz_init_set( &dummy, &(napGetCoeff( p )->z) );
249          term = make_cf( dummy );
250        }
251        else
252        {
253          // assume s==0 or s==1
254          MP_INT num, den;
255          On(SW_RATIONAL);
256          mpz_init_set( &num, &(napGetCoeff( p )->z) );
257          mpz_init_set( &den, &(napGetCoeff( p )->n) );
258          term = make_cf( num, den, ( napGetCoeff( p )->s != 1 ));
259        }
260      }
261    }
262    for ( int i = n; i >0; i-- )
263    {
264      if ( (e = napGetExp( p, i )) != 0 )
265        term *= power( Variable( i ), e );
266    }
267    result += term;
268    p = napNext( p );
269  }
270  return result;
271}
272
273napoly convFactoryPSingTr ( const CanonicalForm & f )
274{
275//    cerr << " f = " << f << endl;
276  int n = napVariables+1;
277  /* ASSERT( level( f ) <= napVariables, "illegal number of variables" ); */
278  int * exp = (int *)omAlloc0(n*sizeof(int));
279  napoly result = NULL;
280  convRecPTr( f, exp, result );
281  omFreeSize((ADDRESS)exp,n*sizeof(int));
282  return result;
283}
284
285static void convRecPTr ( const CanonicalForm & f, int * exp, napoly & result )
286{
287  if (f.isZero())
288    return;
289  if ( ! f.inCoeffDomain() )
290  {
291    int l = f.level();
292    for ( CFIterator i = f; i.hasTerms(); i++ )
293    {
294        exp[l] = i.exp();
295        convRecPTr( i.coeff(), exp, result );
296    }
297    exp[l] = 0;
298  }
299  else
300  {
301    napoly term = napNew();
302    // napNext( term ) = NULL; //already done by napNew
303    for ( int i = 1; i <= napVariables; i++ )
304      napSetExp( term, i , exp[i]);
305    if ( f.isImm() )
306      napGetCoeff( term ) = nacInit( f.intval() );
307    else
308    {
309      number z=(number)omAllocBin(rnumber_bin);
310#if defined(LDEBUG)
311      z->debug=123456;
312#endif
313      z->z = gmp_numerator( f );
314      if ( f.den().isOne() )
315      {
316        z->s = 3;
317      }
318      else
319      {
320        z->n = gmp_denominator( f );
321        if (mpz_cmp_si(&z->n,(long)1)==0)
322        {
323          mpz_clear(&z->n);
324          z->s=3;
325        }
326        else
327        {
328          z->s = 0;
329          nacNormalize(z);
330        }
331      }
332      napGetCoeff( term ) = z;
333    }
334    if (nacIsZero(pGetCoeff(term)))
335    {
336      napDelete(&term);
337    }
338    else
339    {
340      result = napAdd( result, term );
341    }
342  }
343}
344
345CanonicalForm convSingAPFactoryAP ( poly p , const Variable & a)
346{
347  CanonicalForm result = 0;
348  int e, n = pVariables;
349  int off=rPar(currRing);
350
351  On(SW_RATIONAL);
352  while ( p!=NULL)
353  {
354    CanonicalForm term=convSingAFactoryA(((lnumber)pGetCoeff(p))->z,a);
355    for ( int i = 1; i <= n; i++ )
356    {
357      if ( (e = pGetExp( p, i )) != 0 )
358        term *= power( Variable( i + off), e );
359    }
360    result += term;
361    pIter( p );
362  }
363  return result;
364}
365
366static void
367convRecAP_R ( const CanonicalForm & f, int * exp, poly & result, int par_start, int var_start) ;
368
369poly convFactoryAPSingAP_R ( const CanonicalForm & f, int par_start, int var_start )
370{
371  int n = pVariables+1+rPar(currRing);
372  int * exp = (int *)omAlloc0(n*sizeof(int));
373  poly result = NULL;
374  convRecAP_R( f, exp, result,par_start, var_start );
375  omFreeSize((ADDRESS)exp,n*sizeof(int));
376  return result;
377}
378
379poly convFactoryAPSingAP ( const CanonicalForm & f )
380{
381  return convFactoryAPSingAP_R(f,0,rPar(currRing));
382}
383
384static void convRecAP_R ( const CanonicalForm & f, int * exp, poly & result, int par_start, int var_start )
385{
386  if (f.isZero())
387    return;
388  if ( ! f.inCoeffDomain() )
389  {
390    int l = f.level();
391    for ( CFIterator i = f; i.hasTerms(); i++ )
392    {
393      exp[l] = i.exp();
394      convRecAP_R( i.coeff(), exp, result, par_start, var_start );
395    }
396    exp[l] = 0;
397  }
398  else
399  {
400    napoly z=(napoly)convFactoryASingA( f );
401    if (z!=NULL)
402    {
403      poly term = pInit();
404      pNext( term ) = NULL;
405      int i;
406      for ( i = 1; i <= pVariables; i++ )
407        pSetExp( term, i , exp[i+var_start]);
408      //if (rRing_has_Comp(currRing->algring)) p_SetComp(term, 0, currRing->algring); // done by pInit
409      if (par_start==0)
410      {
411        for ( i = 1; i <= var_start; i++ )
412        //z->e[i-1]+=exp[i];
413          napAddExp(z,i,exp[i]);
414      }
415      else
416      {
417        for ( i = par_start+1; i <= var_start+rPar(currRing); i++ )
418        //z->e[i-1]+=exp[i];
419          napAddExp(z,i,exp[i-par_start]);
420      }
421      pGetCoeff(term)=(number)omAlloc0Bin(rnumber_bin);
422      ((lnumber)pGetCoeff(term))->z=z;
423      pSetm( term );
424      result = pAdd( result, term );
425    }
426  }
427}
428
429CanonicalForm convSingAFactoryA ( napoly p , const Variable & a )
430{
431  CanonicalForm result = 0;
432  int e;
433
434  while ( p!=NULL )
435  {
436    CanonicalForm term;
437    /* does only work for finite fields:Z/p */
438    if ( rField_is_Zp_a() )
439    {
440      term = npInt( napGetCoeff( p ) );
441    }
442    else
443    {
444      if ( SR_HDL(napGetCoeff( p )) & SR_INT )
445        term = SR_TO_INT(napGetCoeff( p )) ;
446      else
447      {
448        if ( napGetCoeff( p )->s == 3 )
449        {
450          MP_INT dummy;
451          mpz_init_set( &dummy, &(napGetCoeff( p )->z) );
452          term = make_cf( dummy );
453        }
454        else
455        {
456          // assume s==0 or s==1
457          MP_INT num, den;
458          On(SW_RATIONAL);
459          mpz_init_set( &num, &(napGetCoeff( p )->z) );
460          mpz_init_set( &den, &(napGetCoeff( p )->n) );
461          term = make_cf( num, den, ( napGetCoeff( p )->s != 1 ));
462        }
463      }
464    }
465    if ( (e = napGetExp( p, 1 )) != 0 )
466      term *= power( a , e );
467    result += term;
468    p = napNext( p );
469  }
470  return result;
471}
472
473static number convFactoryNSingAN( const CanonicalForm &f)
474{
475  if ( f.isImm() )
476    return nacInit( f.intval() );
477  else
478  {
479    number z=(number)omAllocBin(rnumber_bin);
480#if defined(LDEBUG)
481    z->debug=123456;
482#endif
483    z->z = gmp_numerator( f );
484    if ( f.den().isOne() )
485    {
486      z->s = 3;
487    }
488    else
489    {
490      z->n = gmp_denominator( f );
491      z->s = 0;
492      nlNormalize(z);
493    }
494    #ifdef LDEBUG
495    nlDBTest(z,__FILE__,__LINE__);
496    #endif
497    return z;
498  }
499}
500
501napoly convFactoryASingA ( const CanonicalForm & f )
502{
503  napoly a=NULL;
504  napoly t;
505  for( CFIterator i=f; i.hasTerms(); i++)
506  {
507    t=napNew();
508    // napNext( t ) = NULL; //already done by napNew
509    napGetCoeff(t)=convFactoryNSingAN( i.coeff() );
510    if (nacIsZero(napGetCoeff(t)))
511    {
512      napDelete(&t);
513    }
514    else
515    {
516      napSetExp(t,1,i.exp());
517      a=napAdd(a,t);
518    }
519  }
520  if (a!=NULL)
521  {
522    if (naMinimalPoly!=NULL)
523    {
524      if (napGetExp(a,1) >= napGetExp(naMinimalPoly,1))
525        a = napRemainder( a, naMinimalPoly);
526    }
527  }
528  return a;
529}
530
531CanonicalForm convSingTrPFactoryP ( poly p )
532{
533  CanonicalForm result = 0;
534  int e, n = pVariables;
535  int offs = rPar(currRing);
536
537  while ( p!=NULL )
538  {
539    nNormalize(pGetCoeff(p));
540    CanonicalForm term=convSingTrFactoryP(((lnumber)pGetCoeff(p))->z);
541
542    if ((((lnumber)pGetCoeff(p))->n!=NULL)
543    && (!errorreported))
544    {
545      WerrorS("conversion error: denominator!= 1");
546    }
547
548    for ( int i = n; i > 0; i-- )
549    {
550      if ( (e = pGetExp( p, i )) != 0 )
551        term = term * power( Variable( i + offs ), e );
552    }
553    result += term;
554    p = pNext( p );
555  }
556  return result;
557}
558
559poly convFactoryPSingTrP ( const CanonicalForm & f )
560{
561  int n = pVariables+1;
562  int * exp = (int*)omAlloc0(n*sizeof(int));
563  poly result = NULL;
564  convRecTrP( f, exp, result , rPar(currRing) );
565  omFreeSize((ADDRESS)exp,n*sizeof(int));
566  return result;
567}
568
569static void
570convRecTrP ( const CanonicalForm & f, int * exp, poly & result , int offs)
571{
572  if (f.isZero())
573    return;
574  if ( f.level() > offs )
575  {
576    int l = f.level();
577    for ( CFIterator i = f; i.hasTerms(); i++ )
578    {
579      exp[l-offs] = i.exp();
580      convRecTrP( i.coeff(), exp, result, offs );
581    }
582    exp[l-offs] = 0;
583  }
584  else
585  {
586    poly term = pInit();
587    pNext( term ) = NULL;
588    for ( int i = 1; i <= pVariables; i++ )
589      pSetExp( term, i ,exp[i]);
590    //if (rRing_has_Comp(currRing)) p_SetComp(term, 0, currRing); // done by pInit
591    pGetCoeff(term)=(number)omAlloc0Bin(rnumber_bin);
592    ((lnumber)pGetCoeff(term))->z=convFactoryPSingTr( f );
593    pSetm( term );
594    result = pAdd( result, term );
595  }
596}
597
598number   nlChineseRemainder(number *x, number *q,int rl)
599// elemenst in the array are x[0..(rl-1)], q[0..(rl-1)]
600{
601#ifdef HAVE_FACTORY
602  CFArray X(rl), Q(rl);
603  int i;
604  for(i=rl-1;i>=0;i--)
605  {
606    X[i]=CanonicalForm(nlInt(x[i]));
607    Q[i]=CanonicalForm(nlInt(q[i]));
608  }
609  CanonicalForm xnew,qnew;
610  chineseRemainder(X,Q,xnew,qnew);
611  number n=convFactoryNSingN(xnew);
612  number p=convFactoryNSingN(qnew);
613  number p2=nlIntDiv(p,nlInit(2));
614  if (nlGreater(n,p2))
615  {
616     number n2=nlSub(n,p);
617     nlDelete(&n,currRing);
618     n=n2;
619  }
620  nlDelete(&p,currRing);
621  nlDelete(&p2,currRing);
622  return n;
623#else
624  WerrorS("not implemented");
625  return nlInit(0);
626#endif
627}
628
629CanonicalForm
630convSingGFFactoryGF( poly p )
631{
632  CanonicalForm result=CanonicalForm(0);
633  int e, n = pVariables;
634
635  while ( p != NULL )
636  {
637    CanonicalForm term;
638    term = make_cf_from_gf( (int)(long)pGetCoeff( p ) );
639    //int * A=(int *)&term;
640    //Print("term=%x, == 0 ?: %d\n",*A, term.isZero());
641    for ( int i = 1; i <= n; i++ )
642    {
643      if ( (e = pGetExp( p, i )) != 0 )
644         term *= power( Variable( i ), e );
645    }
646    result += term;
647    p = pNext( p );
648  }
649  return result;
650}
651
652poly
653convFactoryGFSingGF ( const CanonicalForm & f )
654{
655//    cerr << " f = " << f << endl;
656  int n = pVariables+1;
657  /* ASSERT( level( f ) <= pVariables, "illegal number of variables" ); */
658  int * exp = (int*)omAlloc0(n*sizeof(int));
659  poly result = NULL;
660  convRecGFGF( f, exp, result );
661  omFreeSize((ADDRESS)exp,n*sizeof(int));
662  return result;
663}
664
665static void
666convRecGFGF ( const CanonicalForm & f, int * exp, poly & result )
667{
668  if (f.isZero())
669    return;
670  if ( ! f.inCoeffDomain() )
671  {
672    int l = f.level();
673    for ( CFIterator i = f; i.hasTerms(); i++ )
674    {
675      exp[l] = i.exp();
676      convRecGFGF( i.coeff(), exp, result );
677    }
678    exp[l] = 0;
679  }
680  else
681  {
682    poly term = pInit();
683    pNext( term ) = NULL;
684    for ( int i = 1; i <= pVariables; i++ )
685      pSetExp( term, i, exp[i]);
686    //if (rRing_has_Comp(currRing)) p_SetComp(term, 0, currRing); // done by pInit
687    pGetCoeff( term ) = (number) gf_value (f);
688    pSetm( term );
689    result = pAdd( result, term );
690  }
691}
692
693int convFactoryISingI( const CanonicalForm & f)
694{
695  if (!f.isImm()) WerrorS("int overflow in det");
696  return f.intval();
697}
698#endif
Note: See TracBrowser for help on using the repository browser.