← Back to branch summary

~ubuntu-branches/ubuntu/saucy/gfan/saucy-proposed

« back to all changes in this revision

Viewing changes to singular.cpp

  • Committer: Package Import Robot
  • Author(s): Cédric Boutillier
  • Date: 2013-07-09 10:44:01 UTC
  • mfrom: (2.1.2 experimental)
  • Revision ID: package-import@ubuntu.com-20130709104401-5q66ozz5j5af0dak
Tags: 0.5+dfsg-3
* Upload to unstable.
* modify remove_failing_tests_on_32bits.patch to replace command of
  0009RenderStairCase test with an empty one instead of deleting it.
* remove lintian override about spelling error

Show diffs side-by-side

added added

removed removed

Lines of Context:
singular.cpp
 
1
#include "singularconversion.h"
 
2
 
 
3
#include "polynomialring.h"
 
4
#include "polynomial.h"
 
5
#include "field_rationals.h"
 
6
#include "buchberger.h"
 
7
#include "division.h"
 
8
#include "printer.h"
 
9
#include "log.h"
 
10
#include <iostream>
 
11
 
 
12
 
 
13
void singularBuchberger(PolynomialSet *g, TermOrder const &termOrder)
 
14
{
 
15
  ring R=singularRing(g->getRing());
 
16
 
 
17
  ideal i=singularPolynomialSet(*g);
 
18
 
 
19
  //  ideal j=kStd(i,NULL,testHomog,NULL);
 
20
//  test|=(Sy_bit(OPT_REDSB)|Sy_bit(OPT_REDTAIL)|Sy_bit(OPT_INTSTRATEGY));
 
21
  test|=(Sy_bit(OPT_REDSB)|Sy_bit(OPT_REDTAIL));
 
22
 
 
23
  ideal j=kStd(i,NULL,testHomog,NULL);
 
24
  //idShow(j);
 
25
 
 
26
  idDelete(&i);
 
27
  PolynomialSet ret=fromSingularIdeal(g->getRing(),j);
 
28
  idDelete(&j);
 
29
 
 
30
 
 
31
  freeSingularRing(R);
 
32
 
 
33
  *g=ret;
 
34
  //  return ret;
 
35
}
 
36
 
 
37
/**************************************************************************************************************************/
 
38
 
 
39
#include "groebnerengine.h"
 
40
 
 
41
static ring mySingularRingDegRevLex(PolynomialRing const &r, IntegerVector const &weight)
 
42
{
 
43
  //  FieldRationalsImplementation *q=dynamic_cast<FieldRationalsImplementation>r.getField().implementingObject;
 
44
 
 
45
  ring ret=(ring)omAlloc0(sizeof(sip_sring));
 
46
 
 
47
  if(r.getField().isRationals())
 
48
    {
 
49
      ret->ch=0;
 
50
    }
 
51
  else
 
52
    {
 
53
      assert(0);
 
54
    }
 
55
 
 
56
 
 
57
  ret->N=r.getNumberOfVariables();
 
58
 
 
59
  ret->names=(char**) omAlloc(ret->N*sizeof(char*));
 
60
  for(int i=0;i<ret->N;i++)
 
61
    ret->names[i]=omStrDup(r.getVariableName(i).c_str());
 
62
 
 
63
  ret->order=(int*)omAlloc0(3*sizeof(int));
 
64
  ret->block0=(int*)omAlloc0(3*sizeof(int));
 
65
  ret->block1=(int*)omAlloc0(3*sizeof(int));
 
66
 
 
67
  ret->order[0]=ringorder_wp;//degree revlex
 
68
  // ret->order[0]=ringorder_Wp;//degree lex
 
69
  ret->block0[0]=1;ret->block1[0]=ret->N;
 
70
  ret->order[1]=ringorder_C;
 
71
 
 
72
  ret->wvhdl=(int**)omAlloc0(3*sizeof(int*));
 
73
  ret->wvhdl[0]=(int*)omAlloc(ret->N*sizeof(int));
 
74
  {
 
75
          IntegerVector weight2=weight+(1-weight.min())*IntegerVector::allOnes(weight.size());
 
76
          //debug<<"WEIGHT FOR SINGULAR"<<weight2<<"\n";
 
77
  for(int i=0;i<ret->N;i++)
 
78
    ret->wvhdl[0][i]=weight2[i];//FIXMEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEE
 
79
  }
 
80
  rComplete(ret);
 
81
 // ret->options|=(Sy_bit(OPT_REDSB)|Sy_bit(OPT_REDTAIL)|Sy_bit(OPT_INTSTRATEGY)|Sy_bit(OPT_REDTHROUGH));
 
82
//  cerr<<"ringoptions"<<int(ret->options)<<"\n";
 
83
 
 
84
 // cerr<<"test"<<int(test)<<"\n";
 
85
  rChangeCurrRing(ret);
 
86
 // cerr<<"test"<<int(test)<<"\n";
 
87
 
 
88
  return ret;
 
89
}
 
90
 
 
91
 
 
92
class GroebnerEngineSingular : public GroebnerEngine
 
93
{
 
94
  virtual PolynomialSet groebnerBasis(bool &success, PolynomialSet const &idealGenerators, TermOrder const &termOrder, bool autoreduce)
 
95
  {
 
96
    PolynomialSet ret(idealGenerators.getRing());
 
97
    if(!ret.getRing().getField().isRationals())
 
98
      {
 
99
        success=false;
 
100
        return ret;
 
101
      }
 
102
    if(!idealGenerators.isHomogeneous())
 
103
    {
 
104
        success=false;
 
105
        return ret;
 
106
    }
 
107
    if(dynamic_cast<const WeightReverseLexicographicTermOrder*> (&termOrder))
 
108
      {
 
109
const WeightReverseLexicographicTermOrder *T=dynamic_cast<const WeightReverseLexicographicTermOrder*> (&termOrder);
 
110
 
 
111
ring R=mySingularRingDegRevLex(idealGenerators.getRing(),T->getWeight());
 
112
        ideal i=singularPolynomialSet(idealGenerators);
 
113
        test|=(Sy_bit(OPT_REDSB)|Sy_bit(OPT_REDTAIL)|Sy_bit(OPT_INTSTRATEGY));
 
114
        test|=(Sy_bit(OPT_REDTHROUGH));
 
115
 
 
116
        log2 cerr<<"calling singular\n";
 
117
        //  debug<<"test"<<int(test)<<"\n";
 
118
        ideal j=kStd(i,NULL,testHomog,NULL);
 
119
        log2 cerr<<"returning from singular\n";
 
120
 
 
121
        idDelete(&i);
 
122
        ret=fromSingularIdeal(ret.getRing(),j);
 
123
        idDelete(&j);
 
124
        freeSingularRing(R);
 
125
 
 
126
        ret.markAndScale(termOrder);
 
127
        minimize(&ret);//<---- This is needed for test case 1000 in the suite. Of course it would be nicer if Singular could compute the reduced GB itself.
 
128
        autoReduce_(&ret,termOrder);//------------------------------------------REMOVE AND TRUST SINGULAR
 
129
        if(!isReduced(ret))
 
130
                {
 
131
                        /*
 
132
                         * This test fails on the example
 
133
                         * gfan _tropicaltraverse --symmetry --log2 <hannah.cone
 
134
                         * The monomial d^3*e*f^3*g^8*h^4*i^2 appears as a leading term, but also in a tail.
 
135
                         */
 
136
 
 
137
                debug<<idealGenerators;
 
138
                debug<<ret;
 
139
                  autoReduce_(&ret,termOrder);
 
140
                        debug<<ret;
 
141
                        assert(0);
 
142
                }
 
143
        assert(isMinimal(ret));
 
144
        assert(isReduced(ret));
 
145
        minimize(&ret);
 
146
        if(autoreduce)
 
147
          autoReduce_(&ret,termOrder);
 
148
 
 
149
        success=true;
 
150
        return ret;
 
151
      }
 
152
    success=false;
 
153
    return ret;
 
154
  }
 
155
  virtual PolynomialSet autoReduce(bool &success, PolynomialSet const &idealGenerators)
 
156
  {
 
157
    PolynomialSet ret(idealGenerators.getRing());
 
158
    if(!ret.getRing().getField().isRationals())
 
159
      {
 
160
        success=false;
 
161
        return ret;
 
162
      }
 
163
    if(!idealGenerators.isHomogeneous())
 
164
    {
 
165
        success=false;
 
166
        return ret;
 
167
    }
 
168
 
 
169
    IntegerVector weight=termorderWeight(idealGenerators);
 
170
 
 
171
    ring R=mySingularRingDegRevLex(idealGenerators.getRing(),weight);
 
172
    ideal i=singularPolynomialSet(idealGenerators);
 
173
    test|=(Sy_bit(OPT_REDSB)|Sy_bit(OPT_REDTAIL)|Sy_bit(OPT_INTSTRATEGY));
 
174
 
 
175
    log2 cerr<<"calling singular\n";
 
176
    ideal j=kStd(i,NULL,testHomog,NULL);
 
177
        //    ideal j=kInterRed(i);
 
178
    log2 cerr<<"returning from singular\n";
 
179
 
 
180
    idDelete(&i);
 
181
    ret=fromSingularIdeal(ret.getRing(),j);
 
182
    idDelete(&j);
 
183
    freeSingularRing(R);
 
184
 
 
185
    ret.markAndScale(WeightTermOrder(weight));
 
186
    assert(isReduced(ret));
 
187
 
 
188
    success=true;
 
189
    return ret;
 
190
  }
 
191
  virtual const char* name()
 
192
  {
 
193
    return "singular";
 
194
  }
 
195
};
 
196
 
 
197
static GroebnerEngineSingular groebnerEngineSingular;