Macros | Functions
mpr_inout.h File Reference

Go to the source code of this file.

Macros

#define DEFAULT_DIGITS   30
 
#define MPR_DENSE   1
 
#define MPR_SPARSE   2
 

Functions

BOOLEAN nuUResSolve (leftv res, leftv args)
 solve a multipolynomial system using the u-resultant Input ideal must be 0-dimensional and (currRing->N) == IDELEMS(ideal). More...
 
BOOLEAN nuMPResMat (leftv res, leftv arg1, leftv arg2)
 returns module representing the multipolynomial resultant matrix Arguments 2: ideal i, int k k=0: use sparse resultant matrix of Gelfand, Kapranov and Zelevinsky k=1: use resultant matrix of Macaulay (k=0 is default) More...
 
BOOLEAN nuLagSolve (leftv res, leftv arg1, leftv arg2, leftv arg3)
 find the (complex) roots an univariate polynomial Determines the roots of an univariate polynomial using Laguerres' root-solver. More...
 
BOOLEAN nuVanderSys (leftv res, leftv arg1, leftv arg2, leftv arg3)
 COMPUTE: polynomial p with values given by v at points p1,..,pN derived from p; more precisely: consider p as point in K^n and v as N elements in K, let p1,..,pN be the points in K^n obtained by evaluating all monomials of degree 0,1,...,N at p in lexicographical order, then the procedure computes the polynomial f satisfying f(pi) = v[i] RETURN: polynomial f of degree d. More...
 
BOOLEAN loNewtonP (leftv res, leftv arg1)
 compute Newton Polytopes of input polynomials More...
 
BOOLEAN loSimplex (leftv res, leftv args)
 Implementation of the Simplex Algorithm. More...
 

Macro Definition Documentation

◆ DEFAULT_DIGITS

#define DEFAULT_DIGITS   30

Definition at line 13 of file mpr_inout.h.

◆ MPR_DENSE

#define MPR_DENSE   1

Definition at line 15 of file mpr_inout.h.

◆ MPR_SPARSE

#define MPR_SPARSE   2

Definition at line 16 of file mpr_inout.h.

Function Documentation

◆ loNewtonP()

BOOLEAN loNewtonP ( leftv  res,
leftv  arg1 
)

compute Newton Polytopes of input polynomials

Definition at line 4489 of file ipshell.cc.

4490 {
4491  res->data= (void*)loNewtonPolytope( (ideal)arg1->Data() );
4492  return FALSE;
4493 }
#define FALSE
Definition: auxiliary.h:94
ideal loNewtonPolytope(const ideal id)
Definition: mpr_base.cc:3190
void * data
Definition: subexpr.h:88
void * Data()
Definition: subexpr.cc:1137

◆ loSimplex()

BOOLEAN loSimplex ( leftv  res,
leftv  args 
)

Implementation of the Simplex Algorithm.

For args, see class simplex.

Definition at line 4495 of file ipshell.cc.

4496 {
4497  if ( !(rField_is_long_R(currRing)) )
4498  {
4499  WerrorS("Ground field not implemented!");
4500  return TRUE;
4501  }
4502 
4503  simplex * LP;
4504  matrix m;
4505 
4506  leftv v= args;
4507  if ( v->Typ() != MATRIX_CMD ) // 1: matrix
4508  return TRUE;
4509  else
4510  m= (matrix)(v->CopyD());
4511 
4512  LP = new simplex(MATROWS(m),MATCOLS(m));
4513  LP->mapFromMatrix(m);
4514 
4515  v= v->next;
4516  if ( v->Typ() != INT_CMD ) // 2: m = number of constraints
4517  return TRUE;
4518  else
4519  LP->m= (int)(long)(v->Data());
4520 
4521  v= v->next;
4522  if ( v->Typ() != INT_CMD ) // 3: n = number of variables
4523  return TRUE;
4524  else
4525  LP->n= (int)(long)(v->Data());
4526 
4527  v= v->next;
4528  if ( v->Typ() != INT_CMD ) // 4: m1 = number of <= constraints
4529  return TRUE;
4530  else
4531  LP->m1= (int)(long)(v->Data());
4532 
4533  v= v->next;
4534  if ( v->Typ() != INT_CMD ) // 5: m2 = number of >= constraints
4535  return TRUE;
4536  else
4537  LP->m2= (int)(long)(v->Data());
4538 
4539  v= v->next;
4540  if ( v->Typ() != INT_CMD ) // 6: m3 = number of == constraints
4541  return TRUE;
4542  else
4543  LP->m3= (int)(long)(v->Data());
4544 
4545 #ifdef mprDEBUG_PROT
4546  Print("m (constraints) %d\n",LP->m);
4547  Print("n (columns) %d\n",LP->n);
4548  Print("m1 (<=) %d\n",LP->m1);
4549  Print("m2 (>=) %d\n",LP->m2);
4550  Print("m3 (==) %d\n",LP->m3);
4551 #endif
4552 
4553  LP->compute();
4554 
4555  lists lres= (lists)omAlloc( sizeof(slists) );
4556  lres->Init( 6 );
4557 
4558  lres->m[0].rtyp= MATRIX_CMD; // output matrix
4559  lres->m[0].data=(void*)LP->mapToMatrix(m);
4560 
4561  lres->m[1].rtyp= INT_CMD; // found a solution?
4562  lres->m[1].data=(void*)(long)LP->icase;
4563 
4564  lres->m[2].rtyp= INTVEC_CMD;
4565  lres->m[2].data=(void*)LP->posvToIV();
4566 
4567  lres->m[3].rtyp= INTVEC_CMD;
4568  lres->m[3].data=(void*)LP->zrovToIV();
4569 
4570  lres->m[4].rtyp= INT_CMD;
4571  lres->m[4].data=(void*)(long)LP->m;
4572 
4573  lres->m[5].rtyp= INT_CMD;
4574  lres->m[5].data=(void*)(long)LP->n;
4575 
4576  res->data= (void*)lres;
4577 
4578  return FALSE;
4579 }
sleftv * m
Definition: lists.h:45
Class used for (list of) interpreter objects.
Definition: subexpr.h:82
matrix mapToMatrix(matrix m)
void compute()
#define Print
Definition: emacs.cc:83
Definition: tok.h:95
Definition: lists.h:22
#define FALSE
Definition: auxiliary.h:94
Linear Programming / Linear Optimization using Simplex - Algorithm.
Definition: mpr_numeric.h:194
#define TRUE
Definition: auxiliary.h:98
intvec * zrovToIV()
void WerrorS(const char *s)
Definition: feFopen.cc:24
int Typ()
Definition: subexpr.cc:995
#define omAlloc(size)
Definition: omAllocDecl.h:210
void * data
Definition: subexpr.h:88
ring currRing
Widely used global variable which specifies the current polynomial ring for Singular interpreter and ...
Definition: polys.cc:10
intvec * posvToIV()
BOOLEAN mapFromMatrix(matrix m)
ip_smatrix * matrix
int m
Definition: cfEzgcd.cc:119
leftv next
Definition: subexpr.h:86
INLINE_THIS void Init(int l=0)
const Variable & v
< [in] a sqrfree bivariate poly
Definition: facBivar.h:37
#define MATCOLS(i)
Definition: matpol.h:28
slists * lists
Definition: mpr_numeric.h:146
static BOOLEAN rField_is_long_R(const ring r)
Definition: ring.h:534
int rtyp
Definition: subexpr.h:91
void * Data()
Definition: subexpr.cc:1137
#define MATROWS(i)
Definition: matpol.h:27
int icase
Definition: mpr_numeric.h:201
void * CopyD(int t)
Definition: subexpr.cc:707

◆ nuLagSolve()

BOOLEAN nuLagSolve ( leftv  res,
leftv  arg1,
leftv  arg2,
leftv  arg3 
)

find the (complex) roots an univariate polynomial Determines the roots of an univariate polynomial using Laguerres' root-solver.

Good for polynomials with low and middle degree (<40). Arguments 3: poly arg1 , int arg2 , int arg3 arg2>0: defines precision of fractional part if ground field is Q arg3: number of iterations for approximation of roots (default=2) Returns a list of all (complex) roots of the polynomial arg1

Definition at line 4604 of file ipshell.cc.

4605 {
4606 
4607  poly gls;
4608  gls= (poly)(arg1->Data());
4609  int howclean= (int)(long)arg3->Data();
4610 
4611  if ( !(rField_is_R(currRing) ||
4612  rField_is_Q(currRing) ||
4615  {
4616  WerrorS("Ground field not implemented!");
4617  return TRUE;
4618  }
4619 
4622  {
4623  unsigned long int ii = (unsigned long int)arg2->Data();
4624  setGMPFloatDigits( ii, ii );
4625  }
4626 
4627  if ( gls == NULL || pIsConstant( gls ) )
4628  {
4629  WerrorS("Input polynomial is constant!");
4630  return TRUE;
4631  }
4632 
4633  int ldummy;
4634  int deg= currRing->pLDeg( gls, &ldummy, currRing );
4635  int i,vpos=0;
4636  poly piter;
4637  lists elist;
4638  lists rlist;
4639 
4640  elist= (lists)omAlloc( sizeof(slists) );
4641  elist->Init( 0 );
4642 
4643  if ( rVar(currRing) > 1 )
4644  {
4645  piter= gls;
4646  for ( i= 1; i <= rVar(currRing); i++ )
4647  if ( pGetExp( piter, i ) )
4648  {
4649  vpos= i;
4650  break;
4651  }
4652  while ( piter )
4653  {
4654  for ( i= 1; i <= rVar(currRing); i++ )
4655  if ( (vpos != i) && (pGetExp( piter, i ) != 0) )
4656  {
4657  WerrorS("The input polynomial must be univariate!");
4658  return TRUE;
4659  }
4660  pIter( piter );
4661  }
4662  }
4663 
4664  rootContainer * roots= new rootContainer();
4665  number * pcoeffs= (number *)omAlloc( (deg+1) * sizeof( number ) );
4666  piter= gls;
4667  for ( i= deg; i >= 0; i-- )
4668  {
4669  if ( piter && pTotaldegree(piter) == i )
4670  {
4671  pcoeffs[i]= nCopy( pGetCoeff( piter ) );
4672  //nPrint( pcoeffs[i] );PrintS(" ");
4673  pIter( piter );
4674  }
4675  else
4676  {
4677  pcoeffs[i]= nInit(0);
4678  }
4679  }
4680 
4681 #ifdef mprDEBUG_PROT
4682  for (i=deg; i >= 0; i--)
4683  {
4684  nPrint( pcoeffs[i] );PrintS(" ");
4685  }
4686  PrintLn();
4687 #endif
4688 
4689  roots->fillContainer( pcoeffs, NULL, 1, deg, rootContainer::onepoly, 1 );
4690  roots->solver( howclean );
4691 
4692  int elem= roots->getAnzRoots();
4693  char *dummy;
4694  int j;
4695 
4696  rlist= (lists)omAlloc( sizeof(slists) );
4697  rlist->Init( elem );
4698 
4700  {
4701  for ( j= 0; j < elem; j++ )
4702  {
4703  rlist->m[j].rtyp=NUMBER_CMD;
4704  rlist->m[j].data=(void *)nCopy((number)(roots->getRoot(j)));
4705  //rlist->m[j].data=(void *)(number)(roots->getRoot(j));
4706  }
4707  }
4708  else
4709  {
4710  for ( j= 0; j < elem; j++ )
4711  {
4712  dummy = complexToStr( (*roots)[j], gmp_output_digits, currRing->cf );
4713  rlist->m[j].rtyp=STRING_CMD;
4714  rlist->m[j].data=(void *)dummy;
4715  }
4716  }
4717 
4718  elist->Clean();
4719  //omFreeSize( (ADDRESS) elist, sizeof(slists) );
4720 
4721  // this is (via fillContainer) the same data as in root
4722  //for ( i= deg; i >= 0; i-- ) nDelete( &pcoeffs[i] );
4723  //omFreeSize( (ADDRESS) pcoeffs, (deg+1) * sizeof( number ) );
4724 
4725  delete roots;
4726 
4727  res->rtyp= LIST_CMD;
4728  res->data= (void*)rlist;
4729 
4730  return FALSE;
4731 }
complex root finder for univariate polynomials based on laguers algorithm
Definition: mpr_numeric.h:65
sleftv * m
Definition: lists.h:45
void PrintLn()
Definition: reporter.cc:310
Definition: lists.h:22
#define FALSE
Definition: auxiliary.h:94
static BOOLEAN rField_is_R(const ring r)
Definition: ring.h:510
static short rVar(const ring r)
#define rVar(r) (r->N)
Definition: ring.h:583
#define TRUE
Definition: auxiliary.h:98
bool solver(const int polishmode=PM_NONE)
Definition: mpr_numeric.cc:449
void WerrorS(const char *s)
Definition: feFopen.cc:24
static number & pGetCoeff(poly p)
return an alias to the leading coefficient of p assumes that p != NULL NOTE: not copy ...
Definition: monomials.h:51
#define omAlloc(size)
Definition: omAllocDecl.h:210
#define nPrint(a)
only for debug, over any initalized currRing
Definition: numbers.h:46
void * data
Definition: subexpr.h:88
#define pIter(p)
Definition: monomials.h:44
ring currRing
Widely used global variable which specifies the current polynomial ring for Singular interpreter and ...
Definition: polys.cc:10
#define pGetExp(p, i)
Exponent.
Definition: polys.h:41
int j
Definition: myNF.cc:70
static long pTotaldegree(poly p)
Definition: polys.h:264
#define pIsConstant(p)
like above, except that Comp might be != 0
Definition: polys.h:221
void fillContainer(number *_coeffs, number *_ievpoint, const int _var, const int _tdg, const rootType _rt, const int _anz)
Definition: mpr_numeric.cc:312
int i
Definition: cfEzgcd.cc:123
void PrintS(const char *s)
Definition: reporter.cc:284
static BOOLEAN rField_is_Q(const ring r)
Definition: ring.h:501
gmp_complex * getRoot(const int i)
Definition: mpr_numeric.h:88
static BOOLEAN rField_is_long_C(const ring r)
Definition: ring.h:537
INLINE_THIS void Init(int l=0)
#define NULL
Definition: omList.c:10
slists * lists
Definition: mpr_numeric.h:146
int getAnzRoots()
Definition: mpr_numeric.h:97
static BOOLEAN rField_is_long_R(const ring r)
Definition: ring.h:534
int rtyp
Definition: subexpr.h:91
#define nCopy(n)
Definition: numbers.h:15
void Clean(ring r=currRing)
Definition: lists.h:25
void * Data()
Definition: subexpr.cc:1137
Definition: tok.h:117
char * complexToStr(gmp_complex &c, const unsigned int oprec, const coeffs src)
Definition: mpr_complex.cc:706
size_t gmp_output_digits
Definition: mpr_complex.cc:44
void setGMPFloatDigits(size_t digits, size_t rest)
Set size of mantissa digits - the number of output digits (basis 10) the size of mantissa consists of...
Definition: mpr_complex.cc:62
polyrec * poly
Definition: hilb.h:10
#define nInit(i)
Definition: numbers.h:24

◆ nuMPResMat()

BOOLEAN nuMPResMat ( leftv  res,
leftv  arg1,
leftv  arg2 
)

returns module representing the multipolynomial resultant matrix Arguments 2: ideal i, int k k=0: use sparse resultant matrix of Gelfand, Kapranov and Zelevinsky k=1: use resultant matrix of Macaulay (k=0 is default)

Definition at line 4581 of file ipshell.cc.

4582 {
4583  ideal gls = (ideal)(arg1->Data());
4584  int imtype= (int)(long)arg2->Data();
4585 
4586  uResultant::resMatType mtype= determineMType( imtype );
4587 
4588  // check input ideal ( = polynomial system )
4589  if ( mprIdealCheck( gls, arg1->Name(), mtype, true ) != mprOk )
4590  {
4591  return TRUE;
4592  }
4593 
4594  uResultant *resMat= new uResultant( gls, mtype, false );
4595  if (resMat!=NULL)
4596  {
4597  res->rtyp = MODUL_CMD;
4598  res->data= (void*)resMat->accessResMat()->getMatrix();
4599  if (!errorreported) delete resMat;
4600  }
4601  return errorreported;
4602 }
Base class for solving 0-dim poly systems using u-resultant.
Definition: mpr_base.h:62
resMatrixBase * accessResMat()
Definition: mpr_base.h:78
#define TRUE
Definition: auxiliary.h:98
uResultant::resMatType determineMType(int imtype)
const char * Name()
Definition: subexpr.h:120
Definition: mpr_base.h:98
void * data
Definition: subexpr.h:88
virtual ideal getMatrix()
Definition: mpr_base.h:31
mprState mprIdealCheck(const ideal theIdeal, const char *name, uResultant::resMatType mtype, BOOLEAN rmatrix=false)
short errorreported
Definition: feFopen.cc:23
#define NULL
Definition: omList.c:10
int rtyp
Definition: subexpr.h:91
void * Data()
Definition: subexpr.cc:1137

◆ nuUResSolve()

BOOLEAN nuUResSolve ( leftv  res,
leftv  args 
)

solve a multipolynomial system using the u-resultant Input ideal must be 0-dimensional and (currRing->N) == IDELEMS(ideal).

Resultant method can be MPR_DENSE, which uses Macaulay Resultant (good for dense homogeneous polynoms) or MPR_SPARSE, which uses Sparse Resultant (Gelfand, Kapranov, Zelevinsky). Arguments 4: ideal i, int k, int l, int m k=0: use sparse resultant matrix of Gelfand, Kapranov and Zelevinsky k=1: use resultant matrix of Macaulay (k=0 is default) l>0: defines precision of fractional part if ground field is Q m=0,1,2: number of iterations for approximation of roots (default=2) Returns a list containing the roots of the system.

Definition at line 4834 of file ipshell.cc.

4835 {
4836  leftv v= args;
4837 
4838  ideal gls;
4839  int imtype;
4840  int howclean;
4841 
4842  // get ideal
4843  if ( v->Typ() != IDEAL_CMD )
4844  return TRUE;
4845  else gls= (ideal)(v->Data());
4846  v= v->next;
4847 
4848  // get resultant matrix type to use (0,1)
4849  if ( v->Typ() != INT_CMD )
4850  return TRUE;
4851  else imtype= (int)(long)v->Data();
4852  v= v->next;
4853 
4854  if (imtype==0)
4855  {
4856  ideal test_id=idInit(1,1);
4857  int j;
4858  for(j=IDELEMS(gls)-1;j>=0;j--)
4859  {
4860  if (gls->m[j]!=NULL)
4861  {
4862  test_id->m[0]=gls->m[j];
4863  intvec *dummy_w=id_QHomWeight(test_id, currRing);
4864  if (dummy_w!=NULL)
4865  {
4866  WerrorS("Newton polytope not of expected dimension");
4867  delete dummy_w;
4868  return TRUE;
4869  }
4870  }
4871  }
4872  }
4873 
4874  // get and set precision in digits ( > 0 )
4875  if ( v->Typ() != INT_CMD )
4876  return TRUE;
4877  else if ( !(rField_is_R(currRing) || rField_is_long_R(currRing) || \
4879  {
4880  unsigned long int ii=(unsigned long int)v->Data();
4881  setGMPFloatDigits( ii, ii );
4882  }
4883  v= v->next;
4884 
4885  // get interpolation steps (0,1,2)
4886  if ( v->Typ() != INT_CMD )
4887  return TRUE;
4888  else howclean= (int)(long)v->Data();
4889 
4890  uResultant::resMatType mtype= determineMType( imtype );
4891  int i,count;
4892  lists listofroots= NULL;
4893  number smv= NULL;
4894  BOOLEAN interpolate_det= (mtype==uResultant::denseResMat)?TRUE:FALSE;
4895 
4896  //emptylist= (lists)omAlloc( sizeof(slists) );
4897  //emptylist->Init( 0 );
4898 
4899  //res->rtyp = LIST_CMD;
4900  //res->data= (void *)emptylist;
4901 
4902  // check input ideal ( = polynomial system )
4903  if ( mprIdealCheck( gls, args->Name(), mtype ) != mprOk )
4904  {
4905  return TRUE;
4906  }
4907 
4908  uResultant * ures;
4909  rootContainer ** iproots;
4910  rootContainer ** muiproots;
4911  rootArranger * arranger;
4912 
4913  // main task 1: setup of resultant matrix
4914  ures= new uResultant( gls, mtype );
4915  if ( ures->accessResMat()->initState() != resMatrixBase::ready )
4916  {
4917  WerrorS("Error occurred during matrix setup!");
4918  return TRUE;
4919  }
4920 
4921  // if dense resultant, check if minor nonsingular
4922  if ( mtype == uResultant::denseResMat )
4923  {
4924  smv= ures->accessResMat()->getSubDet();
4925 #ifdef mprDEBUG_PROT
4926  PrintS("// Determinant of submatrix: ");nPrint(smv);PrintLn();
4927 #endif
4928  if ( nIsZero(smv) )
4929  {
4930  WerrorS("Unsuitable input ideal: Minor of resultant matrix is singular!");
4931  return TRUE;
4932  }
4933  }
4934 
4935  // main task 2: Interpolate specialized resultant polynomials
4936  if ( interpolate_det )
4937  iproots= ures->interpolateDenseSP( false, smv );
4938  else
4939  iproots= ures->specializeInU( false, smv );
4940 
4941  // main task 3: Interpolate specialized resultant polynomials
4942  if ( interpolate_det )
4943  muiproots= ures->interpolateDenseSP( true, smv );
4944  else
4945  muiproots= ures->specializeInU( true, smv );
4946 
4947 #ifdef mprDEBUG_PROT
4948  int c= iproots[0]->getAnzElems();
4949  for (i=0; i < c; i++) pWrite(iproots[i]->getPoly());
4950  c= muiproots[0]->getAnzElems();
4951  for (i=0; i < c; i++) pWrite(muiproots[i]->getPoly());
4952 #endif
4953 
4954  // main task 4: Compute roots of specialized polys and match them up
4955  arranger= new rootArranger( iproots, muiproots, howclean );
4956  arranger->solve_all();
4957 
4958  // get list of roots
4959  if ( arranger->success() )
4960  {
4961  arranger->arrange();
4962  listofroots= listOfRoots(arranger, gmp_output_digits );
4963  }
4964  else
4965  {
4966  WerrorS("Solver was unable to find any roots!");
4967  return TRUE;
4968  }
4969 
4970  // free everything
4971  count= iproots[0]->getAnzElems();
4972  for (i=0; i < count; i++) delete iproots[i];
4973  omFreeSize( (ADDRESS) iproots, count * sizeof(rootContainer*) );
4974  count= muiproots[0]->getAnzElems();
4975  for (i=0; i < count; i++) delete muiproots[i];
4976  omFreeSize( (ADDRESS) muiproots, count * sizeof(rootContainer*) );
4977 
4978  delete ures;
4979  delete arranger;
4980  nDelete( &smv );
4981 
4982  res->data= (void *)listofroots;
4983 
4984  //emptylist->Clean();
4985  // omFreeSize( (ADDRESS) emptylist, sizeof(slists) );
4986 
4987  return FALSE;
4988 }
int status int void size_t count
Definition: si_signals.h:59
complex root finder for univariate polynomials based on laguers algorithm
Definition: mpr_numeric.h:65
Class used for (list of) interpreter objects.
Definition: subexpr.h:82
void PrintLn()
Definition: reporter.cc:310
Base class for solving 0-dim poly systems using u-resultant.
Definition: mpr_base.h:62
Definition: tok.h:95
Definition: lists.h:22
#define FALSE
Definition: auxiliary.h:94
static BOOLEAN rField_is_R(const ring r)
Definition: ring.h:510
resMatrixBase * accessResMat()
Definition: mpr_base.h:78
#define omFreeSize(addr, size)
Definition: omAllocDecl.h:260
intvec * id_QHomWeight(ideal id, const ring r)
#define TRUE
Definition: auxiliary.h:98
uResultant::resMatType determineMType(int imtype)
void * ADDRESS
Definition: auxiliary.h:115
void pWrite(poly p)
Definition: polys.h:290
void WerrorS(const char *s)
Definition: feFopen.cc:24
int getAnzElems()
Definition: mpr_numeric.h:95
rootContainer ** specializeInU(BOOLEAN matchUp=false, const number subDetVal=NULL)
Definition: mpr_base.cc:3059
int Typ()
Definition: subexpr.cc:995
const char * Name()
Definition: subexpr.h:120
Definition: mpr_base.h:98
#define nPrint(a)
only for debug, over any initalized currRing
Definition: numbers.h:46
void * data
Definition: subexpr.h:88
ring currRing
Widely used global variable which specifies the current polynomial ring for Singular interpreter and ...
Definition: polys.cc:10
Definition: intvec.h:14
int j
Definition: myNF.cc:70
bool success()
Definition: mpr_numeric.h:162
void arrange()
Definition: mpr_numeric.cc:895
int i
Definition: cfEzgcd.cc:123
void PrintS(const char *s)
Definition: reporter.cc:284
void solve_all()
Definition: mpr_numeric.cc:870
#define IDELEMS(i)
Definition: simpleideals.h:24
mprState mprIdealCheck(const ideal theIdeal, const char *name, uResultant::resMatType mtype, BOOLEAN rmatrix=false)
rootContainer ** interpolateDenseSP(BOOLEAN matchUp=false, const number subDetVal=NULL)
Definition: mpr_base.cc:2921
#define nDelete(n)
Definition: numbers.h:16
leftv next
Definition: subexpr.h:86
static BOOLEAN rField_is_long_C(const ring r)
Definition: ring.h:537
ideal idInit(int idsize, int rank)
initialise an ideal / module
Definition: simpleideals.cc:38
const Variable & v
< [in] a sqrfree bivariate poly
Definition: facBivar.h:37
#define nIsZero(n)
Definition: numbers.h:19
#define NULL
Definition: omList.c:10
static BOOLEAN rField_is_long_R(const ring r)
Definition: ring.h:534
void * Data()
Definition: subexpr.cc:1137
size_t gmp_output_digits
Definition: mpr_complex.cc:44
void setGMPFloatDigits(size_t digits, size_t rest)
Set size of mantissa digits - the number of output digits (basis 10) the size of mantissa consists of...
Definition: mpr_complex.cc:62
virtual IStateType initState() const
Definition: mpr_base.h:41
int BOOLEAN
Definition: auxiliary.h:85
lists listOfRoots(rootArranger *self, const unsigned int oprec)
Definition: ipshell.cc:4991
virtual number getSubDet()
Definition: mpr_base.h:37

◆ nuVanderSys()

BOOLEAN nuVanderSys ( leftv  res,
leftv  arg1,
leftv  arg2,
leftv  arg3 
)

COMPUTE: polynomial p with values given by v at points p1,..,pN derived from p; more precisely: consider p as point in K^n and v as N elements in K, let p1,..,pN be the points in K^n obtained by evaluating all monomials of degree 0,1,...,N at p in lexicographical order, then the procedure computes the polynomial f satisfying f(pi) = v[i] RETURN: polynomial f of degree d.

Definition at line 4733 of file ipshell.cc.

4734 {
4735  int i;
4736  ideal p,w;
4737  p= (ideal)arg1->Data();
4738  w= (ideal)arg2->Data();
4739 
4740  // w[0] = f(p^0)
4741  // w[1] = f(p^1)
4742  // ...
4743  // p can be a vector of numbers (multivariate polynom)
4744  // or one number (univariate polynom)
4745  // tdg = deg(f)
4746 
4747  int n= IDELEMS( p );
4748  int m= IDELEMS( w );
4749  int tdg= (int)(long)arg3->Data();
4750 
4751  res->data= (void*)NULL;
4752 
4753  // check the input
4754  if ( tdg < 1 )
4755  {
4756  WerrorS("Last input parameter must be > 0!");
4757  return TRUE;
4758  }
4759  if ( n != rVar(currRing) )
4760  {
4761  Werror("Size of first input ideal must be equal to %d!",rVar(currRing));
4762  return TRUE;
4763  }
4764  if ( m != (int)pow((double)tdg+1,(double)n) )
4765  {
4766  Werror("Size of second input ideal must be equal to %d!",
4767  (int)pow((double)tdg+1,(double)n));
4768  return TRUE;
4769  }
4770  if ( !(rField_is_Q(currRing) /* ||
4771  rField_is_R() || rField_is_long_R() ||
4772  rField_is_long_C()*/ ) )
4773  {
4774  WerrorS("Ground field not implemented!");
4775  return TRUE;
4776  }
4777 
4778  number tmp;
4779  number *pevpoint= (number *)omAlloc( n * sizeof( number ) );
4780  for ( i= 0; i < n; i++ )
4781  {
4782  pevpoint[i]=nInit(0);
4783  if ( (p->m)[i] )
4784  {
4785  tmp = pGetCoeff( (p->m)[i] );
4786  if ( nIsZero(tmp) || nIsOne(tmp) || nIsMOne(tmp) )
4787  {
4788  omFreeSize( (ADDRESS)pevpoint, n * sizeof( number ) );
4789  WerrorS("Elements of first input ideal must not be equal to -1, 0, 1!");
4790  return TRUE;
4791  }
4792  } else tmp= NULL;
4793  if ( !nIsZero(tmp) )
4794  {
4795  if ( !pIsConstant((p->m)[i]))
4796  {
4797  omFreeSize( (ADDRESS)pevpoint, n * sizeof( number ) );
4798  WerrorS("Elements of first input ideal must be numbers!");
4799  return TRUE;
4800  }
4801  pevpoint[i]= nCopy( tmp );
4802  }
4803  }
4804 
4805  number *wresults= (number *)omAlloc( m * sizeof( number ) );
4806  for ( i= 0; i < m; i++ )
4807  {
4808  wresults[i]= nInit(0);
4809  if ( (w->m)[i] && !nIsZero(pGetCoeff((w->m)[i])) )
4810  {
4811  if ( !pIsConstant((w->m)[i]))
4812  {
4813  omFreeSize( (ADDRESS)pevpoint, n * sizeof( number ) );
4814  omFreeSize( (ADDRESS)wresults, m * sizeof( number ) );
4815  WerrorS("Elements of second input ideal must be numbers!");
4816  return TRUE;
4817  }
4818  wresults[i]= nCopy(pGetCoeff((w->m)[i]));
4819  }
4820  }
4821 
4822  vandermonde vm( m, n, tdg, pevpoint, FALSE );
4823  number *ncpoly= vm.interpolateDense( wresults );
4824  // do not free ncpoly[]!!
4825  poly rpoly= vm.numvec2poly( ncpoly );
4826 
4827  omFreeSize( (ADDRESS)pevpoint, n * sizeof( number ) );
4828  omFreeSize( (ADDRESS)wresults, m * sizeof( number ) );
4829 
4830  res->data= (void*)rpoly;
4831  return FALSE;
4832 }
vandermonde system solver for interpolating polynomials from their values
Definition: mpr_numeric.h:28
#define FALSE
Definition: auxiliary.h:94
return P p
Definition: myNF.cc:203
#define omFreeSize(addr, size)
Definition: omAllocDecl.h:260
static short rVar(const ring r)
#define rVar(r) (r->N)
Definition: ring.h:583
#define TRUE
Definition: auxiliary.h:98
#define nIsOne(n)
Definition: numbers.h:25
void * ADDRESS
Definition: auxiliary.h:115
void WerrorS(const char *s)
Definition: feFopen.cc:24
#define nIsMOne(n)
Definition: numbers.h:26
static number & pGetCoeff(poly p)
return an alias to the leading coefficient of p assumes that p != NULL NOTE: not copy ...
Definition: monomials.h:51
#define omAlloc(size)
Definition: omAllocDecl.h:210
void * data
Definition: subexpr.h:88
ring currRing
Widely used global variable which specifies the current polynomial ring for Singular interpreter and ...
Definition: polys.cc:10
int m
Definition: cfEzgcd.cc:119
#define pIsConstant(p)
like above, except that Comp might be != 0
Definition: polys.h:221
int i
Definition: cfEzgcd.cc:123
static BOOLEAN rField_is_Q(const ring r)
Definition: ring.h:501
#define IDELEMS(i)
Definition: simpleideals.h:24
#define nIsZero(n)
Definition: numbers.h:19
#define NULL
Definition: omList.c:10
const CanonicalForm & w
Definition: facAbsFact.cc:55
#define nCopy(n)
Definition: numbers.h:15
void * Data()
Definition: subexpr.cc:1137
polyrec * poly
Definition: hilb.h:10
#define nInit(i)
Definition: numbers.h:24
Rational pow(const Rational &a, int e)
Definition: GMPrat.cc:418
void Werror(const char *fmt,...)
Definition: reporter.cc:189