source: git/kernel/numeric/mpr_inout.cc @ 9f7665

spielwiese
Last change on this file since 9f7665 was 9f7665, checked in by Oleksandr Motsak <motsak@…>, 10 years ago
Removed HAVE_CONFIG guards fix: fixed the inclusion of configure-generated *config.h's
  • Property mode set to 100644
File size: 5.0 KB
Line 
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4
5
6/*
7* ABSTRACT - multipolynomial resultant
8*/
9
10
11
12
13#include <kernel/mod2.h>
14
15//#ifdef HAVE_MPR
16
17//-> includes
18#include <misc/auxiliary.h>
19#include <omalloc/omalloc.h>
20
21#include <misc/mylimits.h>
22#include <misc/options.h>
23#include <misc/intvec.h>
24
25#include <coeffs/numbers.h>
26#include <coeffs/mpr_global.h>
27
28#include <polys/matpol.h>
29
30
31#include <kernel/febase.h>
32#include <kernel/structs.h>
33#include <kernel/polys.h>
34#include <kernel/ideals.h>
35
36
37#include <math.h>
38
39#include "mpr_base.h"
40#include "mpr_numeric.h"
41#include "mpr_inout.h"
42
43// to get detailed timigs, define MPR_TIMING
44#ifdef MPR_TIMING
45#define TIMING
46#endif
47#include <factory/timing.h>
48TIMING_DEFINE_PRINT(mpr_overall)
49TIMING_DEFINE_PRINT(mpr_check)
50TIMING_DEFINE_PRINT(mpr_constr)
51TIMING_DEFINE_PRINT(mpr_ures)
52TIMING_DEFINE_PRINT(mpr_mures)
53TIMING_DEFINE_PRINT(mpr_arrange)
54TIMING_DEFINE_PRINT(mpr_solver)
55
56#define TIMING_EPR(t,msg) TIMING_END_AND_PRINT(t,msg);TIMING_RESET(t);
57
58//<-
59
60//-> nPrint(number n)
61void nPrint(number n)
62{
63  poly o=pOne();
64  pSetCoeff(o, nCopy(n) );
65  pWrite0( o );
66  pDelete( &o );
67}
68//<-
69
70//------------------------------------------------------------------------------
71
72//-> void mprPrintError( mprState state )
73void mprPrintError( mprState state, const char * name )
74{
75  switch (state)
76  {
77  case mprWrongRType:
78    WerrorS("Unknown resultant matrix type choosen!");
79    break;
80  case mprHasOne:
81    Werror("One element of the ideal %s is constant!",name);
82    break;
83  case mprInfNumOfVars:
84    Werror("Wrong number of elements in given ideal %s, should be %d resp. %d!",
85           name,(currRing->N)+1,(currRing->N));
86    break;
87  case mprNotZeroDim:
88    Werror("The given ideal %s must be 0-dimensional!",name);
89    break;
90  case mprNotHomog:
91    Werror("The given ideal %s has to be homogeneous in the first ring variable!",
92           name);
93    break;
94  case mprNotReduced:
95    Werror("The given ideal %s has to reduced!",name);
96    break;
97  case mprUnSupField:
98    WerrorS("Ground field not implemented!");
99    break;
100  default:
101    break;
102  }
103}
104//<-
105
106//-> mprState mprIdealCheck()
107mprState mprIdealCheck( const ideal theIdeal,
108                        const char * /*name*/,
109                        uResultant::resMatType mtype,
110                        BOOLEAN rmatrix )
111{
112  mprState state = mprOk;
113  // int power;
114  int k;
115
116  int numOfVars= mtype == uResultant::denseResMat?(currRing->N)-1:(currRing->N);
117  if ( rmatrix ) numOfVars++;
118
119  if ( mtype == uResultant::none )
120    state= mprWrongRType;
121
122  if ( IDELEMS(theIdeal) != numOfVars )
123    state= mprInfNumOfVars;
124
125  for ( k= IDELEMS(theIdeal) - 1; (state == mprOk) && (k >= 0); k-- )
126  {
127    poly p = (theIdeal->m)[k];
128    if ( pIsConstant(p) ) state= mprHasOne;
129    else
130    if ( (mtype == uResultant::denseResMat) && !p_IsHomogeneous(p, currRing) )
131      state=mprNotHomog;
132  }
133
134  if ( !(rField_is_R(currRing)||
135         rField_is_Q(currRing)||
136         rField_is_long_R(currRing)||
137         rField_is_long_C(currRing)||
138         (rmatrix && rField_is_Q_a(currRing))) )
139    state= mprUnSupField;
140
141  if ( state != mprOk ) mprPrintError( state, "" /* name */ );
142
143  return state;
144}
145//<-
146
147//-> uResultant::resMatType determineMType( int imtype )
148uResultant::resMatType determineMType( int imtype )
149{
150  switch ( imtype )
151  {
152  case MPR_DENSE:
153    return uResultant::denseResMat;
154  case 0:
155  case MPR_SPARSE:
156    return uResultant::sparseResMat;
157  default:
158    return uResultant::none;
159  }
160}
161//<-
162
163//-> function u_resultant_det
164poly u_resultant_det( ideal gls, int imtype )
165{
166  uResultant::resMatType mtype= determineMType( imtype );
167  poly resdet;
168  poly emptypoly= pInit();
169  number smv= NULL;
170
171  TIMING_START(mpr_overall);
172
173  // check input ideal ( = polynomial system )
174  if ( mprIdealCheck( gls, "", mtype ) != mprOk )
175  {
176    return emptypoly;
177  }
178
179  uResultant *ures;
180
181  // main task 1: setup of resultant matrix
182  TIMING_START(mpr_constr);
183  ures= new uResultant( gls, mtype );
184  TIMING_EPR(mpr_constr,"construction");
185
186  // if dense resultant, check if minor nonsingular
187  if ( mtype == uResultant::denseResMat )
188  {
189    smv= ures->accessResMat()->getSubDet();
190#ifdef mprDEBUG_PROT
191    PrintS("// Determinant of submatrix: ");nPrint(smv); PrintLn();
192#endif
193    if ( nIsZero(smv) )
194    {
195      WerrorS("Unsuitable input ideal: Minor of resultant matrix is singular!");
196      return emptypoly;
197    }
198  }
199
200  // main task 2: Interpolate resultant polynomial
201  TIMING_START(mpr_ures);
202  resdet= ures->interpolateDense( smv );
203  TIMING_EPR(mpr_ures,"ures");
204
205  // free mem
206  delete ures;
207  nDelete( &smv );
208  pDelete( &emptypoly );
209
210  TIMING_EPR(mpr_overall,"overall");
211
212  return ( resdet );
213}
214//<-
215
216//-----------------------------------------------------------------------------
217
218//#endif // HAVE_MPR
219
220// local Variables: ***
221// folded-file: t ***
222// compile-command-1: "make installg" ***
223// compile-command-2: "make install" ***
224// End: ***
225
226// in folding: C-c x
227// leave fold: C-c y
228//   foldmode: F10
Note: See TracBrowser for help on using the repository browser.