source: git/kernel/clapconv.cc @ d0f98e

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