Changeset 0758b5 in git


Ignore:
Timestamp:
Jul 2, 2012, 2:03:47 PM (12 years ago)
Author:
Christian Eder
Branches:
(u'spielwiese', 'fe61d9c35bf7c61f2b6cbf1b56e25e2f08d536cc')
Children:
c78bdede4a03cf59fc377dd13948c2c482c883fa
Parents:
fee33eddc902b0ad109daba1a48c0d83f7329667
git-author:
Christian Eder <ederc@mathematik.uni-kl.de>2012-07-02 14:03:47+02:00
git-committer:
Christian Eder <ederc@mathematik.uni-kl.de>2012-07-05 16:12:50+02:00
Message:
spielwiese compiles with sba
Todo: Add sba() call
Location:
kernel
Files:
3 edited

Legend:

Unmodified
Added
Removed

  • kernel/kInline.h

    rfee33e r0758b5  
    11731173}
    11741174
     1175// dummy function for function pointer strat->rewCrit being usable in all
     1176// possible choices for criteria
     1177KINLINE BOOLEAN arriRewDummy(poly sig, unsigned long not_sevSig, kStrategy strat, int start=0)
     1178{
     1179  return FALSE;
     1180}
     1181
    11751182#endif // defined(KINLINE) || defined(KUTIL_CC)
    11761183#endif // KINLINE_H

  • kernel/kspoly.cc

    rfee33e r0758b5  
    155155#endif
    156156 
     157#if defined(KDEBUG) && defined(TEST_OPT_DEBUG_RED)
     158  if (TEST_OPT_DEBUG)
     159  {
     160    Print(" to: "); PR->wrp(); Print("\n");
     161  }
     162#endif
     163  return ret;
     164}
     165
     166/***************************************************************
     167 *
     168 * Reduces PR with PW
     169 * Assumes PR != NULL, PW != NULL, Lm(PW) divides Lm(PR)
     170 *
     171 ***************************************************************/
     172int ksReducePolySig(LObject* PR,
     173                 TObject* PW,
     174                 long idx,
     175                 poly spNoether,
     176                 number *coef,
     177                 kStrategy strat)
     178{
     179#ifdef KDEBUG
     180  red_count++;
     181#ifdef TEST_OPT_DEBUG_RED
     182  if (TEST_OPT_DEBUG)
     183  {
     184    Print("Red %d:", red_count); PR->wrp(); Print(" with:");
     185    PW->wrp();
     186  }
     187#endif
     188#endif
     189  int ret = 0;
     190  ring tailRing = PR->tailRing;
     191  kTest_L(PR);
     192  kTest_T(PW);
     193
     194  // signature-based stuff:
     195  // checking for sig-safeness first
     196  // NOTE: This has to be done in the current ring
     197  //
     198  /**********************************************
     199   *
     200   * TODO:
     201   * --------------------------------------------
     202   * if strat->incremental
     203   * Since we are subdividing lower index and
     204   * current index reductions it is enough to
     205   * look at the polynomial part of the signature
     206   * for a check. This should speed-up checking
     207   * a lot!
     208   * if !strat->incremental
     209   * We are not subdividing lower and current index
     210   * due to the fact that we are using the induced
     211   * Schreyer order
     212   *
     213   * nevertheless, this different behaviour is
     214   * taken care of by is_sigsafe
     215   * => one reduction procedure can be used for
     216   * both, the incremental and the non-incremental
     217   * attempt!
     218   * --------------------------------------------
     219   *
     220   *********************************************/
     221  //printf("COMPARE IDX: %ld -- %ld\n",idx,strat->currIdx);
     222  if (!PW->is_sigsafe)
     223  {
     224    poly f1 = p_Copy(PR->GetLmCurrRing(),currRing);
     225    poly f2 = PW->GetLmCurrRing();
     226    poly sigMult = pCopy(PW->sig);   // copy signature of reducer
     227    p_ExpVectorSub(f1, f2, currRing); // Calculate the Monomial we must multiply to p2
     228//#if 1
     229#ifdef DEBUGF5
     230    printf("IN KSREDUCEPOLYSIG: \n");
     231    pWrite(pHead(f1));
     232    pWrite(pHead(f2));
     233    pWrite(sigMult);
     234    printf("--------------\n");
     235#endif
     236    sigMult = pp_Mult_qq(f1,sigMult,currRing);
     237//#if 1
     238#ifdef DEBUGF5
     239    printf("------------------- IN KSREDUCEPOLYSIG: --------------------\n");
     240    pWrite(pHead(f1));
     241    pWrite(pHead(f2));
     242    pWrite(sigMult);
     243    pWrite(PR->sig);
     244    printf("--------------\n");
     245#endif
     246    int sigSafe = p_LmCmp(PR->sig,sigMult,currRing);
     247    // now we can delete the copied polynomial data used for checking for
     248    // sig-safeness of the reduction step
     249//#if 1
     250#ifdef DEBUGF5
     251    printf("%d -- %d sig\n",sigSafe,PW->is_sigsafe);
     252
     253#endif
     254    pDelete(&f1);
     255    pDelete(&sigMult);
     256    // go on with the computations only if the signature of p2 is greater than the
     257    // signature of fm*p1
     258    if(sigSafe != 1)
     259    {
     260      PR->is_redundant = TRUE;
     261      return 3;
     262    }
     263    PW->is_sigsafe  = TRUE;
     264  }
     265  PR->is_redundant = FALSE;
     266  poly p1 = PR->GetLmTailRing();   // p2 | p1
     267  poly p2 = PW->GetLmTailRing();   // i.e. will reduce p1 with p2; lm = LT(p1) / LM(p2)
     268  poly t2 = pNext(p2), lm = p1;    // t2 = p2 - LT(p2); really compute P = LC(p2)*p1 - LT(p1)/LM(p2)*p2
     269  assume(p1 != NULL && p2 != NULL);// Attention, we have rings and there LC(p2) and LC(p1) are special
     270  p_CheckPolyRing(p1, tailRing);
     271  p_CheckPolyRing(p2, tailRing);
     272
     273  pAssume1(p2 != NULL && p1 != NULL &&
     274      p_DivisibleBy(p2,  p1, tailRing));
     275
     276  pAssume1(p_GetComp(p1, tailRing) == p_GetComp(p2, tailRing) ||
     277      (p_GetComp(p2, tailRing) == 0 &&
     278       p_MaxComp(pNext(p2),tailRing) == 0));
     279
     280#ifdef HAVE_PLURAL
     281  if (rIsPluralRing(currRing))
     282  {
     283    // for the time being: we know currRing==strat->tailRing
     284    // no exp-bound checking needed
     285    // (only needed if exp-bound(tailring)<exp-b(currRing))
     286    if (PR->bucket!=NULL)  nc_kBucketPolyRed(PR->bucket, p2,coef);
     287    else
     288    {
     289      poly _p = (PR->t_p != NULL ? PR->t_p : PR->p);
     290      assume(_p != NULL);
     291      nc_PolyPolyRed(_p, p2, coef, currRing);
     292      if (PR->t_p!=NULL) PR->t_p=_p; else PR->p=_p;
     293      PR->pLength=0; // usaully not used, GetpLength re-comoutes it if needed
     294    }
     295    return 0;
     296  }
     297#endif
     298
     299  if (t2==NULL)           // Divisor is just one term, therefore it will
     300  {                       // just cancel the leading term
     301    PR->LmDeleteAndIter();
     302    if (coef != NULL) *coef = n_Init(1, tailRing);
     303    return 0;
     304  }
     305
     306  p_ExpVectorSub(lm, p2, tailRing); // Calculate the Monomial we must multiply to p2
     307
     308  if (tailRing != currRing)
     309  {
     310    // check that reduction does not violate exp bound
     311    while (PW->max != NULL && !p_LmExpVectorAddIsOk(lm, PW->max, tailRing))
     312    {
     313      // undo changes of lm
     314      p_ExpVectorAdd(lm, p2, tailRing);
     315      if (strat == NULL) return 2;
     316      if (! kStratChangeTailRing(strat, PR, PW)) return -1;
     317      tailRing = strat->tailRing;
     318      p1 = PR->GetLmTailRing();
     319      p2 = PW->GetLmTailRing();
     320      t2 = pNext(p2);
     321      lm = p1;
     322      p_ExpVectorSub(lm, p2, tailRing);
     323      ret = 1;
     324    }
     325  }
     326
     327  // take care of coef buisness
     328  if (! n_IsOne(pGetCoeff(p2), tailRing))
     329  {
     330    number bn = pGetCoeff(lm);
     331    number an = pGetCoeff(p2);
     332    int ct = ksCheckCoeff(&an, &bn, tailRing->cf);    // Calculate special LC
     333    p_SetCoeff(lm, bn, tailRing);
     334    if ((ct == 0) || (ct == 2))
     335      PR->Tail_Mult_nn(an);
     336    if (coef != NULL) *coef = an;
     337    else n_Delete(&an, tailRing);
     338  }
     339  else
     340  {
     341    if (coef != NULL) *coef = n_Init(1, tailRing);
     342  }
     343
     344
     345  // and finally,
     346  PR->Tail_Minus_mm_Mult_qq(lm, t2, PW->GetpLength() - 1, spNoether);
     347  assume(PW->GetpLength() == pLength(PW->p != NULL ? PW->p : PW->t_p));
     348  PR->LmDeleteAndIter();
     349
     350  // the following is commented out: shrinking
     351#ifdef HAVE_SHIFTBBA_NONEXISTENT
     352  if ( (currRing->isLPring) && (!strat->homog) )
     353  {
     354    // assume? h->p in currRing
     355    PR->GetP();
     356    poly qq = p_Shrink(PR->p, currRing->isLPring, currRing);
     357    PR->Clear(); // does the right things
     358    PR->p = qq;
     359    PR->t_p = NULL;
     360    PR->SetShortExpVector();
     361  }
     362#endif
     363
    157364#if defined(KDEBUG) && defined(TEST_OPT_DEBUG_RED)
    158365  if (TEST_OPT_DEBUG)

  • kernel/kutil.cc

    rfee33e r0758b5  
    17451745  if (strat->interred_flag) return;
    17461746
    1747   int      l,j,compare;
     1747  int      l;
    17481748  poly m1 = NULL,m2 = NULL; // we need the multipliers for the s-polynomial to compute
    17491749              // the corresponding signatures for criteria checks
    … …  
    18451845  // store from which element this pair comes from for further tests
    18461846  Lp.from = strat->sl+1;   
    1847   if(sigCmp==pOrdSgn)
     1847  if(sigCmp==currRing->OrdSgn)
    18481848  {
    18491849    // pSig > sSig
    … …  
    18861886             // generalized prod-crit for lie-type
    18871887             strat->cp++;
    1888              Lp.p = nc_p_Bracket_qq(pCopy(p),strat->S[i]);
     1888             Lp.p = nc_p_Bracket_qq(pCopy(p),strat->S[i], currRing);
    18891889         }
    18901890         else
    … …  
    44704470{
    44714471if (length<0) return 0;
    4472 if (pLmCmp(set[length].sig,p->sig)== pOrdSgn)
     4472if (pLmCmp(set[length].sig,p->sig)== currRing->OrdSgn)
    44734473  return length+1;
    44744474
    … …  
    44804480  if (an >= en-1)
    44814481  {
    4482     if (pLmCmp(set[an].sig,p->sig) == pOrdSgn) return en;
     4482    if (pLmCmp(set[an].sig,p->sig) == currRing->OrdSgn) return en;
    44834483    return an;
    44844484  }
    44854485  i=(an+en) / 2;
    4486   if (pLmCmp(set[i].sig,p->sig) == pOrdSgn) an=i;
     4486  if (pLmCmp(set[i].sig,p->sig) == currRing->OrdSgn) an=i;
    44874487  else                                      en=i;
    44884488  /*aend. fuer lazy == in !=- machen */
    … …  
    49124912}
    49134913
     4914/*
     4915 * SYZYGY CRITERION for signature-based standard basis algorithms
     4916 */
     4917BOOLEAN syzCriterion(poly sig, unsigned long not_sevSig, kStrategy strat)
     4918{
     4919//#if 1
     4920#ifdef DEBUGF5
     4921  Print("syzygy criterion checks:  ");
     4922  pWrite(sig);
     4923#endif
     4924  for (int k=0; k<strat->syzl; k++)
     4925  {
     4926//#if 1
     4927#ifdef DEBUGF5
     4928    Print("checking with: %d --  ",k);
     4929    pWrite(pHead(strat->syz[k]));
     4930#endif
     4931    if (p_LmShortDivisibleBy(strat->syz[k], strat->sevSyz[k], sig, not_sevSig, currRing))
     4932    {
     4933//#if 1
     4934#ifdef DEBUGF5
     4935      printf("DELETE!\n");
     4936#endif
     4937      return TRUE;
     4938    }
     4939  }
     4940  return FALSE;
     4941}
     4942
     4943/*
     4944 * SYZYGY CRITERION for signature-based standard basis algorithms
     4945 */
     4946BOOLEAN syzCriterionInc(poly sig, unsigned long not_sevSig, kStrategy strat)
     4947{
     4948//#if 1
     4949#ifdef DEBUGF5
     4950  Print("syzygy criterion checks:  ");
     4951  pWrite(sig);
     4952#endif
     4953  int comp = p_GetComp(sig, currRing);
     4954  int min, max;
     4955  if (comp<=1)
     4956    return FALSE;
     4957  else
     4958  {
     4959    min = strat->syzIdx[comp-2];
     4960    //printf("SYZIDX %d/%d\n",strat->syzIdx[comp-2],comp-2);
     4961    //printf("SYZIDX %d/%d\n",strat->syzIdx[comp-1],comp-1);
     4962    //printf("SYZIDX %d/%d\n",strat->syzIdx[comp],comp);
     4963    if (comp == strat->currIdx)
     4964    {
     4965      max = strat->syzl;
     4966    }
     4967    else
     4968    {
     4969      max = strat->syzIdx[comp-1];
     4970    }
     4971    for (int k=min; k<max; k++)
     4972    {
     4973#ifdef DEBUGF5
     4974      printf("COMP %d/%d - MIN %d - MAX %d - SYZL %ld\n",comp,strat->currIdx,min,max,strat->syzl);
     4975      Print("checking with: %d --  ",k);
     4976      pWrite(pHead(strat->syz[k]));
     4977#endif
     4978      if (p_LmShortDivisibleBy(strat->syz[k], strat->sevSyz[k], sig, not_sevSig, currRing))
     4979        return TRUE;
     4980    }
     4981    return FALSE;
     4982  }
     4983}
     4984
     4985/*
     4986 * REWRITTEN CRITERION for signature-based standard basis algorithms
     4987 */
     4988BOOLEAN faugereRewCriterion(poly sig, unsigned long not_sevSig, kStrategy strat, int start=0)
     4989{
     4990  //printf("Faugere Rewritten Criterion\n");
     4991//#if 1
     4992#ifdef DEBUGF5
     4993  printf("rewritten criterion checks:  ");
     4994  pWrite(sig);
     4995#endif
     4996  //for(int k = start; k<strat->sl+1; k++)
     4997  for(int k = strat->sl; k>start; k--)
     4998  {
     4999//#if 1
     5000#ifdef DEBUGF5
     5001    Print("checking with:  ");
     5002    pWrite(strat->sig[k]);
     5003    pWrite(pHead(strat->S[k]));
     5004#endif
     5005    if (p_LmShortDivisibleBy(strat->sig[k], strat->sevSig[k], sig, not_sevSig, currRing))
     5006    //if (p_LmEqual(strat->sig[k], sig, currRing))
     5007    {
     5008//#if 1
     5009#ifdef DEBUGF5
     5010      printf("DELETE!\n");
     5011#endif
     5012      return TRUE;
     5013    }
     5014  }
     5015#ifdef DEBUGF5
     5016  Print("ALL ELEMENTS OF S\n----------------------------------------\n");
     5017  for(int kk = 0; kk<strat->sl+1; kk++)
     5018  {
     5019    pWrite(pHead(strat->S[kk]));
     5020  }
     5021  Print("------------------------------\n");
     5022#endif
     5023  return FALSE;
     5024}
     5025
     5026/*
     5027 * REWRITTEN CRITERION for signature-based standard basis algorithms
     5028 ***************************************************************************
     5029 * TODO:This should become the version of Arri/Perry resp. Bjarke/Stillman *
     5030 ***************************************************************************
     5031 */
     5032
     5033// real implementation of arri's rewritten criterion, only called once in
     5034// kstd2.cc, right before starting reduction
     5035// IDEA:  Arri says that it is enough to consider 1 polynomial for each unique
     5036//        signature appearing during the computations. Thus we first of all go
     5037//        through strat->L and delete all other pairs of the same signature,
     5038//        keeping only the one with least possible leading monomial. After this
     5039//        we check if we really need to compute this critical pair at all: There
     5040//        can be elements already in strat->S whose signatures divide the
     5041//        signature of the critical pair in question and whose multiplied
     5042//        leading monomials are smaller than the leading monomial of the
     5043//        critical pair. In this situation we can discard the critical pair
     5044//        completely.
     5045BOOLEAN arriRewCriterion(poly sig, unsigned long not_sevSig, kStrategy strat, int start=0)
     5046{
     5047  //printf("Arri Rewritten Criterion\n");
     5048  while (strat->Ll > 0 && pLmEqual(strat->L[strat->Ll].sig,strat->P.sig))
     5049  {
     5050    // deletes the short spoly
     5051#ifdef HAVE_RINGS
     5052    if (rField_is_Ring(currRing))
     5053      pLmDelete(strat->L[strat->Ll].p);
     5054    else
     5055#endif
     5056      pLmFree(strat->L[strat->Ll].p);
     5057
     5058    // TODO: needs some masking
     5059    // TODO: masking needs to vanish once the signature
     5060    //       sutff is completely implemented
     5061    strat->L[strat->Ll].p = NULL;
     5062    poly m1 = NULL, m2 = NULL;
     5063
     5064    // check that spoly creation is ok
     5065    while (strat->tailRing != currRing &&
     5066          !kCheckSpolyCreation(&(strat->L[strat->Ll]), strat, m1, m2))
     5067    {
     5068      assume(m1 == NULL && m2 == NULL);
     5069      // if not, change to a ring where exponents are at least
     5070      // large enough
     5071      if (!kStratChangeTailRing(strat))
     5072      {
     5073        WerrorS("OVERFLOW...");
     5074        break;
     5075      }
     5076    }
     5077    // create the real one
     5078    ksCreateSpoly(&(strat->L[strat->Ll]), NULL, strat->use_buckets,
     5079                  strat->tailRing, m1, m2, strat->R);
     5080    if (strat->P.GetLmCurrRing() == NULL)
     5081    {
     5082      deleteInL(strat->L,&strat->Ll,strat->Ll,strat);
     5083    }
     5084    if (strat->L[strat->Ll].GetLmCurrRing() == NULL)
     5085    {
     5086      strat->P.Delete();
     5087      strat->P = strat->L[strat->Ll];
     5088      strat->Ll--;
     5089    }
     5090
     5091    if (strat->P.GetLmCurrRing() != NULL && strat->L[strat->Ll].GetLmCurrRing() != NULL)
     5092    {
     5093      if (pLmCmp(strat->P.GetLmCurrRing(),strat->L[strat->Ll].GetLmCurrRing()) == -1)
     5094      {
     5095        deleteInL(strat->L,&strat->Ll,strat->Ll,strat);
     5096      }
     5097      else
     5098      {
     5099        strat->P.Delete();
     5100        strat->P = strat->L[strat->Ll];
     5101        strat->Ll--;
     5102      }
     5103    }
     5104  }
     5105  for (int ii=strat->sl; ii>-1; ii--)
     5106  {
     5107    if (p_LmShortDivisibleBy(strat->sig[ii], strat->sevSig[ii], strat->P.sig, ~strat->P.sevSig, currRing))
     5108    {
     5109      if (!(pLmCmp(ppMult_mm(strat->P.sig,pHead(strat->S[ii])),ppMult_mm(strat->sig[ii],strat->P.GetLmCurrRing())) == 1))
     5110      {
     5111        strat->P.Delete();
     5112        return TRUE;
     5113      }
     5114    }
     5115  }
     5116  return FALSE;
     5117}
     5118
    49145119/***************************************************************
    49155120 *
    … …  
    56045809        LObject h;
    56055810        h.p = pCopy(Q->m[i]);
    5606         if (pOrdSgn==-1)
     5811        if (currRing->OrdSgn==-1)
    56075812        {
    56085813          deleteHC(&h,strat);
    … …  
    56595864      if (h.p!=NULL)
    56605865      {
    5661         if (pOrdSgn==-1)
     5866        if (currRing->OrdSgn==-1)
    56625867        {
    56635868          cancelunit(&h);  /*- tries to cancel a unit -*/
    … …  
    57175922  if ((strat->Ll>=0)
    57185923#ifdef HAVE_RINGS
    5719        && nIsUnit(pGetCoeff(strat->L[strat->Ll].p))
     5924       && n_IsUnit(pGetCoeff(strat->L[strat->Ll].p), currRing->cf)
    57205925#endif
    57215926       && pIsConstant(strat->L[strat->Ll].p))
    … …  
    58786083  strat->S=strat->Shdl->m;
    58796084
     6085  /*- put polys into S -*/
     6086  if (Q!=NULL)
     6087  {
     6088    strat->fromQ=initec(i);
     6089    memset(strat->fromQ,0,i*sizeof(int));
     6090    for (i=0; i<IDELEMS(Q); i++)
     6091    {
     6092      if (Q->m[i]!=NULL)
     6093      {
     6094        LObject h;
     6095        h.p = pCopy(Q->m[i]);
     6096        //if (TEST_OPT_INTSTRATEGY)
     6097        //{
     6098        //  //pContent(h.p);
     6099        //  h.pCleardenom(); // also does a pContent
     6100        //}
     6101        //else
     6102        //{
     6103        //  h.pNorm();
     6104        //}
     6105        if (currRing->OrdSgn==-1)
     6106        {
     6107          deleteHC(&h,strat);
     6108        }
     6109        if (h.p!=NULL)
     6110        {
     6111          strat->initEcart(&h);
     6112          if (strat->sl==-1)
     6113            pos =0;
     6114          else
     6115          {
     6116            pos = posInS(strat,strat->sl,h.p,h.ecart);
     6117          }
     6118          h.sev = pGetShortExpVector(h.p);
     6119          strat->enterS(h,pos,strat, strat->tl+1);
     6120          enterT(h, strat);
     6121          strat->fromQ[pos]=1;
     6122        }
     6123      }
     6124    }
     6125  }
     6126  /*- put polys into S -*/
     6127  for (i=0; i<IDELEMS(F); i++)
     6128  {
     6129    if (F->m[i]!=NULL)
     6130    {
     6131      LObject h;
     6132      h.p = pCopy(F->m[i]);
     6133      if (currRing->OrdSgn==-1)
     6134      {
     6135        deleteHC(&h,strat);
     6136      }
     6137      else
     6138      {
     6139        h.p=redtailBba(h.p,strat->sl,strat);
     6140      }
     6141      if (h.p!=NULL)
     6142      {
     6143        strat->initEcart(&h);
     6144        if (strat->sl==-1)
     6145          pos =0;
     6146        else
     6147          pos = posInS(strat,strat->sl,h.p,h.ecart);
     6148        h.sev = pGetShortExpVector(h.p);
     6149        strat->enterS(h,pos,strat, strat->tl+1);
     6150        enterT(h,strat);
     6151      }
     6152    }
     6153  }
     6154  for (i=0; i<IDELEMS(P); i++)
     6155  {
     6156    if (P->m[i]!=NULL)
     6157    {
     6158      LObject h;
     6159      h.p=pCopy(P->m[i]);
     6160      if (TEST_OPT_INTSTRATEGY)
     6161      {
     6162        h.pCleardenom();
     6163      }
     6164      else
     6165      {
     6166        h.pNorm();
     6167      }
     6168      if(strat->sl>=0)
     6169      {
     6170        if (currRing->OrdSgn==1)
     6171        {
     6172          h.p=redBba(h.p,strat->sl,strat);
     6173          if (h.p!=NULL)
     6174          {
     6175            h.p=redtailBba(h.p,strat->sl,strat);
     6176          }
     6177        }
     6178        else
     6179        {
     6180          h.p=redMora(h.p,strat->sl,strat);
     6181        }
     6182        if(h.p!=NULL)
     6183        {
     6184          strat->initEcart(&h);
     6185          if (TEST_OPT_INTSTRATEGY)
     6186          {
     6187            h.pCleardenom();
     6188          }
     6189          else
     6190          {
     6191            h.is_normalized = 0;
     6192            h.pNorm();
     6193          }
     6194          h.sev = pGetShortExpVector(h.p);
     6195          h.SetpFDeg();
     6196          pos = posInS(strat,strat->sl,h.p,h.ecart);
     6197          enterpairsSpecial(h.p,strat->sl,h.ecart,pos,strat,strat->tl+1);
     6198          strat->enterS(h,pos,strat, strat->tl+1);
     6199          enterT(h,strat);
     6200        }
     6201      }
     6202      else
     6203      {
     6204        h.sev = pGetShortExpVector(h.p);
     6205        strat->initEcart(&h);
     6206        strat->enterS(h,0,strat, strat->tl+1);
     6207        enterT(h,strat);
     6208      }
     6209    }
     6210  }
     6211}
     6212/*2
     6213*construct the set s from F and {P}
     6214*/
     6215
     6216void initSSpecialSba (ideal F, ideal Q, ideal P,kStrategy strat)
     6217{
     6218  int   i,pos;
     6219
     6220  if (Q!=NULL) i=((IDELEMS(Q)+(setmaxTinc-1))/setmaxTinc)*setmaxTinc;
     6221  else i=setmaxT;
     6222  i=((i+IDELEMS(F)+IDELEMS(P)+15)/16)*16;
     6223  strat->fromS=initec(i);
     6224  strat->sevS=initsevS(i);
     6225  strat->sevSig=initsevS(i);
     6226  strat->S_2_R=initS_2_R(i);
     6227  strat->fromQ=NULL;
     6228  strat->Shdl=idInit(i,F->rank);
     6229  strat->S=strat->Shdl->m;
     6230  strat->sig=(poly *)omAlloc0(i*sizeof(poly));
    58806231  /*- put polys into S -*/
    58816232  if (Q!=NULL)
    … …  
    70817432}
    70827433
     7434void initSbaPos (kStrategy strat)
     7435{
     7436  if (currRing->OrdSgn==1)
     7437  {
     7438    if (strat->honey)
     7439    {
     7440      strat->posInL = posInL15;
     7441      // ok -- here is the deal: from my experiments for Singular-2-0
     7442      // I conclude that that posInT_EcartpLength is the best of
     7443      // posInT15, posInT_EcartFDegpLength, posInT_FDegLength, posInT_pLength
     7444      // see the table at the end of this file
     7445      if (TEST_OPT_OLDSTD)
     7446        strat->posInT = posInT15;
     7447      else
     7448        strat->posInT = posInT_EcartpLength;
     7449    }
     7450    else if (currRing->pLexOrder && !TEST_OPT_INTSTRATEGY)
     7451    {
     7452      strat->posInL = posInL11;
     7453      strat->posInT = posInT11;
     7454    }
     7455    else if (TEST_OPT_INTSTRATEGY)
     7456    {
     7457      strat->posInL = posInL11;
     7458      strat->posInT = posInT11;
     7459    }
     7460    else
     7461    {
     7462      strat->posInL = posInL0;
     7463      strat->posInT = posInT0;
     7464    }
     7465    //if (strat->minim>0) strat->posInL =posInLSpecial;
     7466    if (strat->homog)
     7467    {
     7468       strat->posInL = posInL110;
     7469       strat->posInT = posInT110;
     7470    }
     7471  }
     7472  else
     7473  {
     7474    if (strat->homog)
     7475    {
     7476      strat->posInL = posInL11;
     7477      strat->posInT = posInT11;
     7478    }
     7479    else
     7480    {
     7481      if ((currRing->order[0]==ringorder_c)
     7482      ||(currRing->order[0]==ringorder_C))
     7483      {
     7484        strat->posInL = posInL17_c;
     7485        strat->posInT = posInT17_c;
     7486      }
     7487      else
     7488      {
     7489        strat->posInL = posInL17;
     7490        strat->posInT = posInT17;
     7491      }
     7492    }
     7493  }
     7494  if (strat->minim>0) strat->posInL =posInLSpecial;
     7495  // for further tests only
     7496  if ((BTEST1(11)) || (BTEST1(12)))
     7497    strat->posInL = posInL11;
     7498  else if ((BTEST1(13)) || (BTEST1(14)))
     7499    strat->posInL = posInL13;
     7500  else if ((BTEST1(15)) || (BTEST1(16)))
     7501    strat->posInL = posInL15;
     7502  else if ((BTEST1(17)) || (BTEST1(18)))
     7503    strat->posInL = posInL17;
     7504  if (BTEST1(11))
     7505    strat->posInT = posInT11;
     7506  else if (BTEST1(13))
     7507    strat->posInT = posInT13;
     7508  else if (BTEST1(15))
     7509    strat->posInT = posInT15;
     7510  else if ((BTEST1(17)))
     7511    strat->posInT = posInT17;
     7512  else if ((BTEST1(19)))
     7513    strat->posInT = posInT19;
     7514  else if (BTEST1(12) || BTEST1(14) || BTEST1(16) || BTEST1(18))
     7515    strat->posInT = posInT1;
     7516#ifdef HAVE_RINGS
     7517  if (rField_is_Ring(currRing))
     7518  {
     7519    strat->posInL = posInL11;
     7520    strat->posInT = posInT11;
     7521  }
     7522#endif
     7523  strat->posInLDependsOnLength = FALSE;
     7524  strat->posInLSba  = posInLSig;
     7525  //strat->posInL     = posInLSig;
     7526  strat->posInL     = posInLF5C;
     7527  //strat->posInT     = posInTSig;
     7528}
     7529
     7530void initSbaBuchMora (ideal F,ideal Q,kStrategy strat)
     7531{
     7532  strat->interpt = BTEST1(OPT_INTERRUPT);
     7533  strat->kHEdge=NULL;
     7534  if (currRing->OrdSgn==1) strat->kHEdgeFound=FALSE;
     7535  /*- creating temp data structures------------------- -*/
     7536  strat->cp = 0;
     7537  strat->c3 = 0;
     7538  strat->tail = pInit();
     7539  /*- set s -*/
     7540  strat->sl = -1;
     7541  /*- set ps -*/
     7542  strat->syzl = -1;
     7543  /*- set L -*/
     7544  strat->Lmax = ((IDELEMS(F)+setmaxLinc-1)/setmaxLinc)*setmaxLinc;
     7545  strat->Ll = -1;
     7546  strat->L = initL(((IDELEMS(F)+setmaxLinc-1)/setmaxLinc)*setmaxLinc);
     7547  /*- set B -*/
     7548  strat->Bmax = setmaxL;
     7549  strat->Bl = -1;
     7550  strat->B = initL();
     7551  /*- set T -*/
     7552  strat->tl = -1;
     7553  strat->tmax = setmaxT;
     7554  strat->T = initT();
     7555  strat->R = initR();
     7556  strat->sevT = initsevT();
     7557  /*- init local data struct.---------------------------------------- -*/
     7558  strat->P.ecart=0;
     7559  strat->P.length=0;
     7560  if (currRing->OrdSgn==-1)
     7561  {
     7562    if (strat->kHEdge!=NULL) pSetComp(strat->kHEdge, strat->ak);
     7563    if (strat->kNoether!=NULL) pSetComp(strat->kNoetherTail(), strat->ak);
     7564  }
     7565  if(TEST_OPT_SB_1)
     7566  {
     7567    int i;
     7568    ideal P=idInit(IDELEMS(F)-strat->newIdeal,F->rank);
     7569    for (i=strat->newIdeal;i<IDELEMS(F);i++)
     7570    {
     7571      P->m[i-strat->newIdeal] = F->m[i];
     7572      F->m[i] = NULL;
     7573    }
     7574    initSSpecialSba(F,Q,P,strat);
     7575    for (i=strat->newIdeal;i<IDELEMS(F);i++)
     7576    {
     7577      F->m[i] = P->m[i-strat->newIdeal];
     7578      P->m[i-strat->newIdeal] = NULL;
     7579    }
     7580    idDelete(&P);
     7581  }
     7582  else
     7583  {
     7584    /*Shdl=*/initSLSba(F, Q,strat); /*sets also S, ecartS, fromQ */
     7585    // /*Shdl=*/initS(F, Q,strat); /*sets also S, ecartS, fromQ */
     7586  }
     7587  strat->fromT = FALSE;
     7588  strat->noTailReduction = !TEST_OPT_REDTAIL;
     7589  if (!TEST_OPT_SB_1)
     7590  {
     7591    updateS(TRUE,strat);
     7592  }
     7593  if (strat->fromQ!=NULL) omFreeSize(strat->fromQ,IDELEMS(strat->Shdl)*sizeof(int));
     7594  strat->fromQ=NULL;
     7595}
     7596
     7597void exitSba (kStrategy strat)
     7598{
     7599  /*- release temp data -*/
     7600  cleanT(strat);
     7601  omFreeSize(strat->T,(strat->tmax)*sizeof(TObject));
     7602  omFreeSize(strat->R,(strat->tmax)*sizeof(TObject*));
     7603  omFreeSize(strat->sevT, (strat->tmax)*sizeof(unsigned long));
     7604  omFreeSize(strat->ecartS,IDELEMS(strat->Shdl)*sizeof(int));
     7605  omFreeSize(strat->fromS,IDELEMS(strat->Shdl)*sizeof(int));
     7606  omFreeSize((ADDRESS)strat->sevS,IDELEMS(strat->Shdl)*sizeof(unsigned long));
     7607  omFreeSize((ADDRESS)strat->sevSig,IDELEMS(strat->Shdl)*sizeof(unsigned long));
     7608  omFreeSize((ADDRESS)strat->syz,(strat->syzmax)*sizeof(poly));
     7609  omFreeSize((ADDRESS)strat->sevSyz,(strat->syzmax)*sizeof(unsigned long));
     7610  if (strat->incremental)
     7611  {
     7612    omFreeSize(strat->syzIdx,(strat->syzidxmax)*sizeof(int));
     7613  }
     7614  omFreeSize(strat->S_2_R,IDELEMS(strat->Shdl)*sizeof(int));
     7615  /*- set L: should be empty -*/
     7616  omFreeSize(strat->L,(strat->Lmax)*sizeof(LObject));
     7617  /*- set B: should be empty -*/
     7618  omFreeSize(strat->B,(strat->Bmax)*sizeof(LObject));
     7619  /*- set sig: no need for the signatures anymore -*/
     7620  omFreeSize(strat->sig,IDELEMS(strat->Shdl)*sizeof(poly));
     7621  pLmDelete(&strat->tail);
     7622  strat->syzComp=0;
     7623}
     7624
    70837625/*2
    70847626* in the case of a standardbase of a module over a qring:
    … …  
    75478089
    75488090  kStratChangeTailRing(strat, NULL, NULL, e);
     8091}
     8092
     8093ring sbaRing (kStrategy strat, const ring r, BOOLEAN complete, int sgn)
     8094{
     8095  int n = rBlocks(r); // Including trailing zero!
     8096  // if incremental => use (C,monomial order from r)
     8097  if (strat->incremental)
     8098  {
     8099    if (r->order[0] == ringorder_C || r->order[0] == ringorder_c)
     8100    {
     8101      return r;
     8102    }
     8103      ring res = rCopy0(r, FALSE, TRUE);
     8104      for (int i=1; i<n-1; i++)
     8105      {
     8106        res->order[i] = res->order[i-1];
     8107        res->block0[i] = res->block0[i-1];
     8108        res->block1[i] = res->block1[i-1];
     8109        res->wvhdl[i] = res->wvhdl[i-1];
     8110      }
     8111
     8112    // new 1st block
     8113    res->order[0]   = ringorder_C; // Prefix
     8114    res->block0[0]  = 1;
     8115    res->block1[0]  = res->N;
     8116    //res->wvhdl[j]   = NULL;
     8117    // res->order [j] = 0; // The End!
     8118    rComplete(res, 1);
     8119#ifdef HAVE_PLURAL
     8120    if (rIsPluralRing(r))
     8121    {
     8122      if ( nc_rComplete(r, res, false) ) // no qideal!
     8123      {
     8124#ifndef NDEBUG
     8125        WarnS("error in nc_rComplete");
     8126#endif
     8127        // cleanup?
     8128
     8129        //      rDelete(res);
     8130        //      return r;
     8131
     8132        // just go on..
     8133      }
     8134    }
     8135#endif
     8136  strat->tailRing = res;
     8137  return (res);
     8138  }
     8139  // not incremental => use Schreyer order
     8140  // this is done by a trick when initializing the signatures
     8141  // in initSLSba():
     8142  // Instead of using the signature 1e_i for F->m[i], we start
     8143  // with the signature LM(F->m[i])e_i for F->m[i]. Doing this we get a
     8144  // Schreyer order w.r.t. the underlying monomial order.
     8145  // => we do not need to change the underlying polynomial ring at all!
     8146
     8147
     8148  /*
     8149  else
     8150  {
     8151    ring res = rCopy0(r, FALSE, FALSE);
     8152    // Create 2 more blocks for prefix/suffix:
     8153    res->order=(int *)omAlloc0((n+2)*sizeof(int)); // 0  ..  n+1
     8154    res->block0=(int *)omAlloc0((n+2)*sizeof(int));
     8155    res->block1=(int *)omAlloc0((n+2)*sizeof(int));
     8156    int ** wvhdl =(int **)omAlloc0((n+2)*sizeof(int**));
     8157
     8158    // Encapsulate all existing blocks between induced Schreyer ordering markers: prefix and suffix!
     8159    // Note that prefix and suffix have the same ringorder marker and only differ in block[] parameters!
     8160
     8161    // new 1st block
     8162    int j = 0;
     8163    res->order[j] = ringorder_IS; // Prefix
     8164    res->block0[j] = res->block1[j] = 0;
     8165    // wvhdl[j] = NULL;
     8166    j++;
     8167
     8168    for(int i = 0; (i < n) && (r->order[i] != 0); i++, j++) // i = [0 .. n-1] <- non-zero old blocks
     8169    {
     8170      res->order [j] = r->order [i];
     8171      res->block0[j] = r->block0[i];
     8172      res->block1[j] = r->block1[i];
     8173
     8174      if (r->wvhdl[i] != NULL)
     8175      {
     8176        wvhdl[j] = (int*) omMemDup(r->wvhdl[i]);
     8177      } // else wvhdl[j] = NULL;
     8178    }
     8179
     8180    // new last block
     8181    res->order [j] = ringorder_IS; // Suffix
     8182    res->block0[j] = res->block1[j] = sgn; // Sign of v[o]: 1 for C, -1 for c
     8183    // wvhdl[j] = NULL;
     8184    j++;
     8185
     8186    // res->order [j] = 0; // The End!
     8187    res->wvhdl = wvhdl;
     8188
     8189    // j == the last zero block now!
     8190    assume(j == (n+1));
     8191    assume(res->order[0]==ringorder_IS);
     8192    assume(res->order[j-1]==ringorder_IS);
     8193    assume(res->order[j]==0);
     8194
     8195    if (complete)
     8196    {
     8197      rComplete(res, 1);
     8198
     8199#ifdef HAVE_PLURAL
     8200      if (rIsPluralRing(r))
     8201      {
     8202        if ( nc_rComplete(r, res, false) ) // no qideal!
     8203        {
     8204        }
     8205      }
     8206      assume(rIsPluralRing(r) == rIsPluralRing(res));
     8207#endif
     8208
     8209
     8210#ifdef HAVE_PLURAL
     8211      ring old_ring = r;
     8212
     8213#endif
     8214
     8215      if (r->qideal!=NULL)
     8216      {
     8217        res->qideal= idrCopyR_NoSort(r->qideal, r, res);
     8218
     8219        assume(idRankFreeModule(res->qideal, res) == 0);
     8220
     8221#ifdef HAVE_PLURAL
     8222        if( rIsPluralRing(res) )
     8223          if( nc_SetupQuotient(res, r, true) )
     8224          {
     8225            //          WarnS("error in nc_SetupQuotient"); // cleanup?      rDelete(res);       return r;  // just go on...?
     8226          }
     8227
     8228#endif
     8229        assume(idRankFreeModule(res->qideal, res) == 0);
     8230      }
     8231
     8232#ifdef HAVE_PLURAL
     8233      assume((res->qideal==NULL) == (old_ring->qideal==NULL));
     8234      assume(rIsPluralRing(res) == rIsPluralRing(old_ring));
     8235      assume(rIsSCA(res) == rIsSCA(old_ring));
     8236      assume(ncRingType(res) == ncRingType(old_ring));
     8237#endif
     8238    }
     8239    strat->tailRing = res;
     8240    return res;
     8241  }
     8242  */
    75498243}
    75508244
Note: See TracChangeset for help on using the changeset viewer.