Changeset 2986f2 in git


Ignore:
Timestamp:
Jan 18, 2019, 11:50:24 PM (5 years ago)
Author:
Viktor Levandovskyy <levandov@…>
Branches:
(u'spielwiese', 'fe61d9c35bf7c61f2b6cbf1b56e25e2f08d536cc')
Children:
4acf3262d5fd8ab635bdc4a814e5951853e47272
Parents:
078a4a90f17cc237bc78fddcc5ef815a8ed3bc37
Message:
letterplace implementation of the factorization
File:
1 edited

Legend:

Unmodified
Added
Removed
  • Singular/LIB/ncfactor.lib

    r078a4a r2986f2  
    11///////////////////////////////////////////////////////////
    2 version="version ncfactor.lib 4.1.1.0 Jan_2018 "; // $Id$
     2version = "$Id$";
    33category="Noncommutative";
    44info="
    55LIBRARY: ncfactor.lib  Tools for factorization in some noncommutative algebras
    6 AUTHORS: Albert Heinle,     aheinle@uwaterloo.ca
    7 @*       Viktor Levandovskyy,     levandov@math.rwth-aachen.de
    8 
    9 OVERVIEW: New methods for factorization on polynomials
    10   are implemented for three types of algebras, all generated by
    11    2*n, n in NN, generators (n'th Weyl, n'th shift and graded polynomials in n'th q-Weyl algebras)
    12    over a field K.
     6AUTHORS: Albert Heinle,     aheinle at uwaterloo.ca
     7@*       Viktor Levandovskyy,     levandov at math.rwth-aachen.de
     8
     9OVERVIEW: In this library, new methods for factorization on polynomials
     10@* are implemented for several types of algebras, namely
     11@* - finitely presented (and also free) associative algebras (Letterplace subsystem)
     12@* - G-algebras, including (q)-Weyl and (q)-shift algebras in 2n variables
     13@* The determination of the best algorithm available for users input is done
     14@* automatically in the procedure ncfactor().
     15
    1316@* More detailled description of the algorithms and related publications can be found at
    14    @url{https://cs.uwaterloo.ca/\~aheinle/}.
     17@url{https://cs.uwaterloo.ca/\~aheinle/}.
    1518
    1619PROCEDURES:
    17   ncfactor(h);               Factorization in any G-algebra.
     20  ncfactor(h);               Factorization in any finitely presented algebra (incl. G-algebras)
    1821  facWeyl(h);                Factorization in the n'th Weyl algebra
    1922  facFirstWeyl(h);           Factorization in the first Weyl algebra
     
    3740LIB "solve.lib"; //for solve
    3841LIB "poly.lib"; //for content
     42LIB "fpadim.lib"; //for letterplace
    3943
    4044proc tst_ncfactor()
     
    128132  }
    129133  print("Successful.");
    130   print("Testing Min");
    131   if (!test_Min())
    132   {
    133     ERROR("Min is not working properly.");
    134   }
    135   print("Successful.");
    136134  print("Testing isListEqual");
    137135  if (!test_isListEqual())
     
    410408  }
    411409  print("Successful.");
    412   print("Testing delete_duplicates_noteval_and_sort");
     410  print("Testing getMaxDegreeVecLetterPlace");
     411  if(!test_getMaxDegreeVecLetterPlace())
     412  {
     413    ERROR("getMaxDegreeVecLetterPlace is not working properly.");
     414  }
     415  print("Successful.");
     416  print("Testing wordsWithMaxAppearance");
     417  if(!test_wordsWithMaxAppearance())
     418  {
     419    ERROR("wordsWithMaxAppearance is not working properly.");
     420  }
     421  print("Successful.");
     422  print("Testing monsSmallerThanLetterPlace");
     423  if(!test_monsSmallerThanLetterPlace())
     424  {
     425    ERROR("monsSmallerThanLetterPlace is not working properly.");
     426  }
     427  print("Successful.");
     428  print("Testing ncfactor_letterplace_poly_s");
     429  if(!test_ncfactor_letterplace_poly_s())
     430  {
     431    ERROR("ncfactor_letterplace_poly_s is not working properly.");
     432  }
     433  print("Successful.");
     434  print("Testing ncfactor_without_opt_letterplace");
     435  if(!test_ncfactor_without_opt_letterplace())
     436  {
     437    ERROR("ncfactor_without_opt_letterplace is not working properly.");
     438  }
     439  print("Successful.");
    413440  print("All tests ran successfully.");
    414441}
     
    20412068  return(result);
    20422069}//testing increment_intvec
    2043 
    2044 
    2045 //I know that Singular has this function in the package
    2046 //qhmoduli.lib. But I don't see the reason why I should clutter my
    2047 //namespace with functions to study Moduli Spaces of
    2048 //Semi-Quasihomogeneous Singularities (in case somebody wonders)
    2049 static proc Min(def lst)
    2050 "USAGE: Min(lst); list is an iterable element (list/ideal/intvec)
    2051         containing ordered elements.
    2052 PURPOSE: returns the minimal element in lst
    2053 ASSUME: lst contains at least one element
    2054 "{//proc Min
    2055   def minElem = lst[1];
    2056   int i;
    2057   for (i = 2; i<=size(lst);i++)
    2058   {
    2059     if (lst[i]<minElem)
    2060     {
    2061       minElem = lst[i];
    2062     }
    2063   }
    2064   return (minElem);
    2065 }//proc Min
    2066 
    2067 
    2068 static proc test_Min()
    2069 {//Testing Min
    2070   int result = 1;
    2071   //Test 1: One element list
    2072   int expected= 1;
    2073   int obtained = Min(list(1));
    2074   if (expected!=obtained)
    2075   {
    2076     print("Test 1 for Min failed.");
    2077     print("Expected:\n");
    2078     print(expected);
    2079     print("obtained:\n");
    2080     print(obtained);
    2081     result = 0;
    2082   }
    2083   //Test 2: Two element list, first one is Min
    2084   expected = 5;
    2085   obtained = Min(list(5,8));
    2086   if (expected!=obtained)
    2087   {
    2088     print("Test 2 for Min failed.");
    2089     print("Expected:\n");
    2090     print(expected);
    2091     print("obtained:\n");
    2092     print(obtained);
    2093     result = 0;
    2094   }
    2095   //Test 3: Two element list, second one is Min
    2096   expected = 5;
    2097   obtained = Min(list(8,5));
    2098   if (expected!=obtained)
    2099   {
    2100     print("Test 3 for Min failed.");
    2101     print("Expected:\n");
    2102     print(expected);
    2103     print("obtained:\n");
    2104     print(obtained);
    2105     result = 0;
    2106   }
    2107   //Test 4: A lot of elements, Min randomly in between.
    2108   expected = -11;
    2109   obtained = Min(list(5,8,7,2,8,0,-11,0,25));
    2110   if (expected!=obtained)
    2111   {
    2112     print("Test 4 for Min failed.");
    2113     print("Expected:\n");
    2114     print(expected);
    2115     print("obtained:\n");
    2116     print(obtained);
    2117     result = 0;
    2118   }
    2119   return(result);
    2120 }//Testing Min
    2121 
    21222070
    21232071static proc isListEqual(list l1, list l2)
     
    33403288"
    33413289{//isInCommutativeSubRing
     3290  int i; int j; int k;
     3291  intvec tempIntVec;
     3292  int previouslyEncounteredVar = -1;
     3293  if ( attrib(basering, "isLetterplaceRing") >= 1)
     3294  {//separate department: letterplace case
     3295    for (i = 1; i<=size(h); i++)
     3296    {//iterating over the monomials of h
     3297      tempIntVec = lp2iv(h[i]);
     3298      if (size(tempIntVec) > 1)
     3299      {//nontrivial monomial
     3300    j = tempIntVec[1];
     3301    if (previouslyEncounteredVar != -1 and j != previouslyEncounteredVar)
     3302    {//In this case, the previously encountered var is not equal to this one
     3303      return(0);
     3304    }//In this case, the previously encountered var is not equal to this one
     3305    previouslyEncounteredVar = j;
     3306    for (k = 2; k <= size(tempIntVec); k++)
     3307    {//checking if each entry in the intvec is the same
     3308      if (tempIntVec[k]!=j)
     3309      {//more than one variable in monomial
     3310        return(0);
     3311      }//more than one variable in monomial
     3312    }//checking if each entry in the intvec is the same
     3313      }//nontrivial monomial
     3314      else
     3315      {i++; continue;}
     3316    }//iterating over the monomials of h
     3317    return(1);
     3318  }//separate department: letterplace case
    33423319  list tempRingList = ringlist(basering);
    33433320  if (size(tempRingList)<=4)
     
    33463323  }//In this case, the given ring was commutative
    33473324  list appearing_variables;
    3348   int i; int j;
    33493325  intvec degreeIntVec;
    33503326  for (i = 1; i<=nvars(basering);i++)
     
    34263402  {
    34273403    print("Test 4 for isInCommutativeSubRing failed.");
     3404    print("Expected:\n");
     3405    print(expected);
     3406    print("obtained:\n");
     3407    print(obtained);
     3408    result = 0;
     3409  }
     3410  return(result);
     3411  //Test 5: Letterplace ring, negative example
     3412  ring r2 = 0,(x,y,z),dp;
     3413  int d =4; // degree bound
     3414  def R2 = makeLetterplaceRing(d);
     3415  setring R2;
     3416  expected = 0;
     3417  obtained = isInCommutativeSubRing(x + y);
     3418  if (expected !=obtained)
     3419  {
     3420    print("Test 5 for isInCommutativeSubRing failed.");
     3421    print("Expected:\n");
     3422    print(expected);
     3423    print("obtained:\n");
     3424    print(obtained);
     3425    result = 0;
     3426  }
     3427  //Test 6: Letterplace ring, positive example
     3428  expected = 1;
     3429  obtained = isInCommutativeSubRing(x + x*x);
     3430  if (expected !=obtained)
     3431  {
     3432    print("Test 6 for isInCommutativeSubRing failed.");
     3433    print("Expected:\n");
     3434    print(expected);
     3435    print("obtained:\n");
     3436    print(obtained);
     3437    result = 0;
     3438  }
     3439  //Test 7: Letterplace ring, previous bug
     3440  expected = 0;
     3441  obtained = isInCommutativeSubRing(x*x - y*y);
     3442  if (expected !=obtained)
     3443  {
     3444    print("Test 7 for isInCommutativeSubRing failed.");
    34283445    print("Expected:\n");
    34293446    print(expected);
     
    34573474  if (deg(h)<=0)
    34583475  {
    3459     dbprint(p,dbprintWhitespace + "h is a constant. Returning
    3460  immediately");
     3476    dbprint(p,dbprintWhitespace + "h is a constant. Returning immediately");
    34613477    return(list(list(h)));
    34623478  }
    34633479  def r = basering;
    34643480  list factorizeOutput;
    3465   if (size(ringlist(basering))<=4)
     3481  if (size(ringlist(basering))<=4 && ( !(attrib(r, "isLetterplaceRing") >= 1) ) )
    34663482  {//commutative ring case
    3467     dbprint(p,dbprintWhitespace + "We are in a commutative
    3468  ring. Factoring h using factorize.");
     3483    dbprint(p,dbprintWhitespace + "We are in a commutative ring. Factoring h using factorize.");
    34693484    factorizeOutput = factorize(h);
    34703485  }//commutative ring case
    34713486  else
    34723487  {//commutative subring case;
    3473     dbprint(p,dbprintWhitespace + "We are in a commutative
    3474  subring. Generating commutative ring..");
    3475     def rList = ringlist(basering);
    3476     rList = delete(rList,5);
    3477     rList = delete(rList,5);
    3478     def tempRing = ring(rList);
    3479     setring(tempRing);
    3480     poly h = imap(r,h);
    3481     dbprint(p, dbprintWhitespace+"Factoring h in commutative ring.");
    3482     list tempResult = factorize(h);
    3483     setring(r);
    3484     factorizeOutput = imap(tempRing,tempResult);
     3488    dbprint(p,dbprintWhitespace + "We are in a commutative subring. Generating commutative ring.");
     3489    if ( !(attrib(r, "isLetterplaceRing") >= 1) )
     3490    {//G-algebra case
     3491      def rList = ringlist(basering);
     3492      rList = delete(rList,5);
     3493      rList = delete(rList,5);
     3494      def tempRing = ring(rList);
     3495      setring(tempRing);
     3496      poly h = imap(r,h);
     3497      dbprint(p, dbprintWhitespace+"Factoring h in commutative ring.");
     3498      list tempResult = factorize(h);
     3499      setring(r);
     3500      factorizeOutput = imap(tempRing,tempResult);
     3501    }//G-algebra case
     3502    else
     3503    {//Letterplace case
     3504      int theRightVar = lp2iv(h[1])[1];
     3505      list commFactInIv = list();
     3506      list tempMonomList = list();
     3507      def rList = ringlist(basering);
     3508      rList[2] = list("@x");
     3509      rList[3] = list(list("dp", intvec(1)),
     3510              list("C", 0));
     3511      def tempRing = ring(rList);
     3512      setring(tempRing);
     3513      ideal tempIdeal;
     3514      for (i=1; i<=nvars(r); i++)
     3515      {tempIdeal[i] = var(1);}
     3516      map fromR = r, tempIdeal;
     3517      poly h = fromR(h);
     3518      dbprint(p, dbprintWhitespace+"Factoring h in commutative ring.");
     3519      list tempResult = factorize(h);
     3520      dbprint(p, dbprintWhitespace+"Done. The result is:");
     3521      dbprint(p, tempResult);
     3522      dbprint(p, dbprintWhitespace+"Now translated into intvec representation:");
     3523      for (i = 1; i <= size(tempResult[1]); i++)
     3524      {//translate all factorizations into lists of intvecs and coeffs
     3525    commFactInIv[i] = list();
     3526    for (j = 1; j <= size(tempResult[1][i]); j++)
     3527    {//iterate over each monomial
     3528      commFactInIv[i][j] = list(leadcoef(tempResult[1][i][j]));
     3529      if (deg(tempResult[1][i][j]) != 0)
     3530      {//a power of @x
     3531        commFactInIv[i][j] = commFactInIv[i][j] + list(theRightVar:deg(tempResult[1][i][j]));
     3532      }//a power of @x
     3533    }//iterate over each monomial
     3534      }//translate all factorizations into lists of intvecs and coeffs
     3535      dbprint(p, commFactInIv);
     3536      setring(r);
     3537      factorizeOutput = imap(tempRing,tempResult);
     3538      list commFactInIv = imap(tempRing, commFactInIv);
     3539      poly tempEntry;
     3540      for (i = 1; i<=size(commFactInIv); i++)
     3541      {
     3542    tempEntry = 0;
     3543    for (j=1; j<=size(commFactInIv[i]); j++)
     3544    {
     3545      if (size(commFactInIv[i][j]) == 1)
     3546      {//simply a number
     3547        tempEntry = tempEntry + commFactInIv[i][j][1];
     3548      }//simply a number
     3549      else
     3550      {
     3551        tempEntry = tempEntry + commFactInIv[i][j][1]*iv2lp(commFactInIv[i][j][2]);
     3552      }
     3553    }
     3554    factorizeOutput[1][i] = tempEntry;
     3555      }
     3556    }//Letterplace case
    34853557  }//commutative subring case;
    3486   dbprint(p,dbprintWhitespace + "Done commutatively factorizing.
    3487  The result is:");
     3558  dbprint(p,dbprintWhitespace + "Done commutatively factorizing. The result isssssss:");
    34883559  dbprint(p,factorizeOutput);
    34893560  dbprint(p,dbprintWhitespace+"Computing all permutations of this factorization.");
     
    35263597  }
    35273598  kill R;
    3528   //Test 2: Non-commutative ring, constant
     3599  //Test 2: G-Algebra, constant
    35293600  ring R = 0,(x,d),dp;
    35303601  def r = nc_algebra(1,1);
     
    35423613  }
    35433614  kill r; kill R;
    3544   //Test 3: Non-commutative ring, definitely commutative subring, irreducible
     3615  //Test 3: G-Algebra, definitely commutative subring, irreducible
    35453616  ring R = 0,(x1,x2,d1,d2),dp;
    35463617  def r = Weyl();
     
    35583629  }
    35593630  kill r; kill R;
    3560   //Test 4: Non-commutative ring, definitely commutative subring,
     3631  //Test 4: G-Algebra, definitely commutative subring,
    35613632  //reducible
    35623633  ring R = 0,(x1,x2,d1,d2),dp;
     
    35793650  }
    35803651  kill r; kill R;
     3652  //Test 5: Letterplace ring, constant
     3653  ring r = 0,(x,y,z),dp;
     3654  int d =4; // degree bound
     3655  def R = makeLetterplaceRing(d);
     3656  setring R;
     3657  list expected = list(list(poly(3445)));
     3658  list obtained = factor_commutative(3445);
     3659  if (!isListEqual(expected,obtained))
     3660  {
     3661    print("Test 5 for factor_commutative failed.");
     3662    print("Expected:\n");
     3663    print(expected);
     3664    print("obtained:\n");
     3665    print(obtained);
     3666    result = 0;
     3667  }
     3668  //Test 6: Letterplace ring, definitely commutative subring, irreducible
     3669  expected = list(list(poly(1), x*x+ 2*x + 2));
     3670  obtained = factor_commutative(x*x+ 2*x + 2);
     3671  if (!isListEqual(expected,obtained))
     3672  {
     3673    print("Test 6 for factor_commutative failed.");
     3674    print("Expected:\n");
     3675    print(expected);
     3676    print("obtained:\n");
     3677    print(obtained);
     3678    result = 0;
     3679  }
     3680  //Test 7: Letterplace ring, definitely commutative subring,
     3681  //reducible
     3682  expected = list(list(poly(1), x*x+ 2*x + 2, x*x- 2*x - 2),
     3683            list(poly(1), x*x- 2*x - 2, x*x+ 2*x + 2));
     3684  obtained = factor_commutative((x*x+ 2*x + 2)*(x*x- 2*x - 2));
     3685  if (!isListEqual(expected,obtained))
     3686  {
     3687    print("Test 7 for factor_commutative failed.");
     3688    print("Expected:\n");
     3689    print(expected);
     3690    print("obtained:\n");
     3691    print(obtained);
     3692    result = 0;
     3693  }
     3694  //Test 8: Letterplace ring, using last variable
     3695  expected = list(list(poly(1), z*z+ 2*z + 2, z*z - 2*z - 2),
     3696            list(poly(1), z*z - 2*z - 2, z*z+ 2*z + 2));
     3697  obtained = factor_commutative((z*z+ 2*z + 2)*(z*z - 2*z - 2));
     3698  if (!isListEqual(expected,obtained))
     3699  {
     3700    print("Test 8 for factor_commutative failed.");
     3701    print("Expected:\n");
     3702    print(expected);
     3703    print("obtained:\n");
     3704    print(obtained);
     3705    result = 0;
     3706  }
    35813707  return(result);
    35823708}//testing factor_commutative
     
    67306856
    67316857  result = delete_duplicates_noteval_and_sort(result);
    6732   //print(M);
    67336858  if (size(result) == 0)
    67346859  {//no factorization found
     
    68947019  }
    68957020  hZeroinR = tempHList;
    6896   //hBetweenDegrees = reverse(hBetweenDegrees);
    68977021  dbprint(p,dbprintWhitespace+" Done!");
    68987022  dbprint(p,dbprintWhitespace+" Moving everything into the ring K[theta]");
     
    82448368  list tempList = list();
    82458369  //First, we are going to deal with our most hated guys: The Coefficients.
    8246   //
    82478370  dbprint(p,dbprintWhitespace + "We get all combinations for the coefficient of the
    82488371maximal homogeneous part");
     
    97469869  rlist[3][2][2]=intvec(0);
    97479870  rlist = insert(rlist,ideal(0),3);
    9748   //The following lines of code will become obsolete once bug #753 in
    9749   //Singular is fixed. Until then, this measure is unfortunately
    9750   //necessary.
    9751   //begin
    9752   intvec perms = 0:nvars(@r);
    9753   for (i =1; i<=nvars(@r);i++)
    9754   {//finding the right position of each input variable
    9755     for (j = 1; j<=size(inpList); j++)
    9756     {
    9757       if (var(i)==inpList[j])
    9758       {//found the correct spot
    9759     perms[i] = j;
    9760     break;
    9761       }//found the correct spot
    9762     }
    9763   }//finding the right position of each input variable
    97649871  def @r2 = ring(rlist);
    97659872  setring(@r2);
    97669873  def @r3 = Weyl();
    97679874  setring(@r3);
    9768   ideal mappingIdeal;
    9769   for (i = 1; i<=size(perms); i++)
    9770   {//generating the mapping ideal
    9771     if (perms[i]>0)
    9772     {//only variables that are in the input list are permitted
    9773       mappingIdeal[i] = var(perms[i]);
    9774     }//only variables that are in the input list are permitted
    9775   }//generating the mapping ideal
    9776   map toR3 = @r, mappingIdeal;
    9777   poly h = toR3(h);
     9875  poly h = imap(@r,h);
    97789876  list result = facWeyl(h);
    97799877  setring(@r);
    9780   ideal mappingIdeal;
    9781   for (i = 1; i<=size(perms);i++)
    9782   {//finding the correct permutation
    9783     for (j = 1; j<=size(perms); j++)
    9784     {//the spot needs to coincide
    9785       if (perms[j]==i)
    9786       {//found it
    9787     mappingIdeal[i] = var(j);
    9788     break;
    9789       }//found it
    9790     }//the spot needs to coincide
    9791   }//finding the correct permutation
    9792   map toR = @r3, mappingIdeal;
    9793   list result = toR(result);
    9794   //end
    9795   //once bug is fixed, please uncomment the following lines and
    9796   //discard the above
    9797   /* def @r2 = ring(rlist); */
    9798   /* setring(@r2); */
    9799   /* def @r3 = Weyl(); */
    9800   /* setring(@r3); */
    9801   /* poly h = imap(@r,h); */
    9802   /* list result = facWeyl(h); */
    9803   /* setring(@r); */
    9804   /* list result = imap(@r3,result); */
     9878  list result = imap(@r3,result);
    98059879  return(result);
    98069880}//proc facSubWeyl
     
    1047410548 - We have n parameters q_1,..., q_n given.
    1047510549
    10476 SEE ALSO: homogfacFirstQWeyl, homogfacFirstQWeyl_all
     10550SEE ALSO: homogfacNthWeyl, homogfacFirstQWeyl, homogfacFirstQWeyl_all
    1047710551"
    1047810552{//proc homogfacNthQWeyl_all
     
    1100211076  for (i = 1; i<=voice;i++)
    1100311077  {dbprintWhitespace = dbprintWhitespace + " ";}
    11004   dbprint(p,dbprintWhitespace + "Checking if the field fulfills
     11078  dbprint(p,dbprintWhitespace + "Checking if the field fulfills\
    1100511079 everything we assume.");
    1100611080  if (nvars(basering)<1)
     
    1101011084    if (minpoly !=0)
    1101111085    {
    11012       ERROR("factorize does not support fields of type (p^n,a)");
     11086      ERROR("factorize does not support fields of type (p^n,a) yet");
    1101311087      return list(list());
    1101411088    }
    1101511089  }
    11016   dbprint(p,dbprintWhitespace + "Everything seems to be alright with
     11090  dbprint(p,dbprintWhitespace + "Everything seems to be alright with\
    1101711091 the ground field.");
    1101811092  if (deg(h)<=0)
     
    1102411098  dbprint(p,dbprintWhitespace+"Checking if a more improved algorithm\
    1102511099 is available other than the naive ansatz-method.");
    11026   dbprint(p,dbprintWhitespace+"1. Checking if h in commutative
    11027  (sub)-ring");
     11100  dbprint(p,dbprintWhitespace+"1. Checking if h in commutative (sub)-ring");
    1102811101  if (isInCommutativeSubRing(h))
    1102911102  {//h is in a commutative subring
    1103011103    return(factor_commutative(h));
    1103111104  }//h is in a commutative subring
    11032   dbprint(p,dbprintWhitespace+"h was not in a commutative
     11105  dbprint(p,dbprintWhitespace+"h was not in a commutative\
    1103311106 (sub)-ring.");
    11034   dbprint(p,dbprintWhitespace+"2. Checking if h is in Weyl-Algebra
    11035  over QQ");
     11107  dbprint(p,dbprintWhitespace+"2. Checking if h was in a letterplace ring");
     11108  if ( attrib(basering, "isLetterplaceRing") >= 1 )
     11109  {
     11110    dbprint(p,dbprintWhitespace+
     11111            "We are indeed in a letterplace ring. Forwarding the computation " +
     11112            "to ncfactor-without-opt-letterplace");
     11113    return(ncfactor_without_opt_letterplace(h));
     11114  }
     11115  dbprint(p,dbprintWhitespace+"2. Checking if h is in Weyl-Algebra over QQ");
    1103611116  if (ncfactor_isWeyl()&&npars(basering)==0&&char(basering)==0)
    1103711117  {
    11038     dbprint(p,dbprintWhitespace + "We are indeed in a Weyl
    11039  algebra. Forwarding computation to facWeyl");
     11118    dbprint(p,dbprintWhitespace +
     11119            "We are indeed in a Weyl algebra. Forwarding computation to facWeyl");
    1104011120    return(facWeyl(h));
    1104111121  }
    1104211122  dbprint(p,dbprintWhitespace+"Our ring is not a Weyl algebra over QQ.");
    11043   dbprint(p,dbprintWhitespace+"3. Checking if we are in a q-Weyl
    11044  algebra over QQ.");
     11123  dbprint(p,dbprintWhitespace+"3. Checking if we are in a q-Weyl algebra over QQ.");
    1104511124  if (ncfactor_isQWeyl()&&char(basering)==0)
    1104611125  {
    11047     dbprint(p,dbprintWhitespace + "We are indeed in a q-Weyl
     11126    dbprint(p,dbprintWhitespace + "We are indeed in a q-Weyl\
    1104811127 algebra over QQ. Checking if homogeneous.");
    1104911128    if(homogwithorderNthWeyl(h))
    1105011129    {
    11051       dbprint(p,dbprintWhitespace+"h was homogeneous. Forwarding
     11130      dbprint(p,dbprintWhitespace+"h was homogeneous. Forwarding\
    1105211131 computation to homogfacqweyl.");
    1105311132      return(homogfacNthQWeyl_all(h));
     
    1105511134  }
    1105611135  dbprint(p,dbprintWhitespace+"Our ring is not a q-Weyl algebra.");
    11057   dbprint(p,dbprintWhitespace+"4. Checking if h is in a subalgebra
     11136  dbprint(p,dbprintWhitespace+"4. Checking if h is in a subalgebra\
    1105811137 that resembles the Weyl algebra");
    1105911138  if(isInWeylSubAlg(h)&&npars(basering)==0&&char(basering)==0)
    1106011139  {
    11061     dbprint(p,dbprintWhitespace+"We are indeed in a subalgebra that is
     11140    dbprint(p,dbprintWhitespace+"We are indeed in a subalgebra that is\
    1106211141 isomorphic to a Weyl algebra. Forwarding to facSubWeyl.");
    1106311142    return(facSubWeyl(h));
    1106411143  }
    11065   dbprint(p,dbprintWhitespace+"No optimized algorithm available. Going
     11144  dbprint(p,dbprintWhitespace+"No optimized algorithm available. Going \
    1106611145 for the ansatz method without optimization.");
    1106711146  return(ncfactor_without_opt(h));
     
    1128611365  {
    1128711366    print("Test 10 for ncfactor failed.");
     11367    print("Expected:\n");
     11368    print(expected);
     11369    print("obtained:\n");
     11370    print(obtained);
     11371    result = 0;
     11372  }
     11373  kill r;
     11374  //Test 11: Generic Letterplace example
     11375  ring r = 0,(x,y,z),dp;
     11376  int d =4; // degree bound
     11377  def R = makeLetterplaceRing(d);
     11378  setring R;
     11379  poly f1 = 6*x*y*x + 9*x;
     11380  list obtained = ncfactor(f1);
     11381  list expected = sortFactorizations(
     11382    list(
     11383      list(number(3), x, 2*y*x + 3),
     11384      list(number(3), 2*x*y + 3, x)
     11385      ));
     11386  if (!isListEqual(expected,obtained))
     11387  {
     11388    print("Test 11 for ncfactor failed.");
    1128811389    print("Expected:\n");
    1128911390    print(expected);
     
    1198412085}//testing factorize_nc_s
    1198512086
     12087
     12088static proc getMaxDegreeVecLetterPlace(poly h)
     12089"USAGE: getMaxDegreeVecLetterPlace(h); h is a polynomial in a
     12090        Letterplace ring.
     12091RETURN: intvec
     12092PURPOSE: Returns for each variable in the ring the maximal
     12093         degree in which it appears in h.
     12094ASSUMING: The basering is
     12095"
     12096{//proc getMaxDegreeVecLetterPlace
     12097  int lv = attrib(basering, "isLetterplaceRing"); // nvars of orig ring
     12098  if (h == 0) { return (0:lv); }
     12099  intvec tempIntVec1 = 0:lv;
     12100  intvec maxDegrees = 0:lv;
     12101  intvec tempIntVec2;
     12102  int i; int j;
     12103  for (i = 1; i<=size(h); i++)
     12104  {//going through each monomial in h and finding the degree in each var
     12105    tempIntVec2 = lp2iv(h[i]);
     12106    if (tempIntVec2 == 0)
     12107    {
     12108      i++; continue;
     12109    }
     12110    tempIntVec1 = 0:lv;
     12111    for (j=1; j<=size(tempIntVec2); j++)
     12112    {//filling in the number of occurrences
     12113      tempIntVec1[tempIntVec2[j]] = tempIntVec1[tempIntVec2[j]] + 1;
     12114    }//filling in the number of occurrences
     12115    for (j = 1; j<=lv; j++)
     12116    {//putting the max into maxDegrees
     12117      maxDegrees[j] = max(maxDegrees[j], tempIntVec1[j]);
     12118    }//putting the max into maxDegrees
     12119  }//going through each monomial in h and finding the degree in each var
     12120  return(maxDegrees);
     12121}//proc getMaxDegreeVecLetterPlace
     12122
     12123static proc test_getMaxDegreeVecLetterPlace()
     12124{//testing getMaxDegreeVecLetterPlace
     12125  int result = 1;
     12126  ring r = 0,(x,y,z),dp;
     12127  int d =4; // degree bound
     12128  def R = makeLetterplaceRing(d);
     12129  setring R;
     12130  //Test 1: 0
     12131  intvec expected = 0:3;
     12132  intvec obtained = getMaxDegreeVecLetterPlace(0);
     12133  if (expected!=obtained)
     12134  {
     12135    print("Test 1 for getMaxDegreeVecLetterPlace failed.");
     12136    print("Expected:\n");
     12137    print(expected);
     12138    print("obtained:\n");
     12139    print(obtained);
     12140    result = 0;
     12141  }
     12142  //Test 2: Constant neq 0
     12143  expected = 0:3;
     12144  obtained = getMaxDegreeVecLetterPlace(3);
     12145  if (expected!=obtained)
     12146  {
     12147    print("Test 2 for getMaxDegreeVecLetterPlace failed.");
     12148    print("Expected:\n");
     12149    print(expected);
     12150    print("obtained:\n");
     12151    print(obtained);
     12152    result = 0;
     12153  }
     12154  //Test 3: univariate, first variable
     12155  expected = 2, 0, 0;
     12156  obtained = getMaxDegreeVecLetterPlace(5*x*x + x + 1);
     12157  if (expected!=obtained)
     12158  {
     12159    print("Test 3 for getMaxDegreeVecLetterPlace failed.");
     12160    print("Expected:\n");
     12161    print(expected);
     12162    print("obtained:\n");
     12163    print(obtained);
     12164    result = 0;
     12165  }
     12166  //Test 4: univariate, another variable
     12167  expected = 0, 0, 3;
     12168  obtained = getMaxDegreeVecLetterPlace(5*z*z*z + z*z + 1);
     12169  if (expected!=obtained)
     12170  {
     12171    print("Test 4 for getMaxDegreeVecLetterPlace failed.");
     12172    print("Expected:\n");
     12173    print(expected);
     12174    print("obtained:\n");
     12175    print(obtained);
     12176    result = 0;
     12177  }
     12178  //Test 5: random polynomial
     12179  expected = 3, 1, 2;
     12180  obtained = getMaxDegreeVecLetterPlace(2*x*y*x*x + 4*y*z*z);
     12181  if (expected!=obtained)
     12182  {
     12183    print("Test 5 for getMaxDegreeVecLetterPlace failed.");
     12184    print("Expected:\n");
     12185    print(expected);
     12186    print("obtained:\n");
     12187    print(obtained);
     12188    result = 0;
     12189  }
     12190  return(result);
     12191}//testing getMaxDegreeVecLetterPlace
     12192
     12193static proc wordsWithMaxAppearance(intvec maxDegInCoordinates, int upToDeg)
     12194"USAGE: wordsWithMaxAppearance(maxDegInCoordinates);
     12195maxDegreeInCoordinates is an intvec.
     12196RETURN: list
     12197PURPOSE: maxDegInCoordinates is a vector representing
     12198the maximum times a variable is allowed to appear in a monomial in
     12199letterplace representation. This function computes all possible words
     12200given this restriction.
     12201ASSUME: maxDegInCoordinates only has non-negative entries
     12202"{//wordsWithMaxAppearance
     12203  int p=printlevel-voice+2;//for dbprint
     12204  int i; int j;
     12205  string dbprintWhitespace = "";
     12206  for (i = 1; i<=voice;i++)
     12207  {dbprintWhitespace = dbprintWhitespace + " ";}
     12208  list result;
     12209  intvec tempMaxDegs;
     12210  list recResult;
     12211  for (i = 1; i<=size(maxDegInCoordinates); i++)
     12212  {//For each coordinate not equal to zero, do a recursive call and concatenate
     12213    if (maxDegInCoordinates[i] <= 0)
     12214    {//Done with coordinate i
     12215      i++; continue;
     12216    }//Done with coordinate i
     12217    tempMaxDegs = maxDegInCoordinates;
     12218    tempMaxDegs[i] = tempMaxDegs[i] - 1;
     12219    recResult = wordsWithMaxAppearance(tempMaxDegs, upToDeg);
     12220    if (size(recResult) == 0 && upToDeg>=1)
     12221    {//Single entry just needs to be there
     12222      result = result + list(intvec(i));
     12223    }//Single entry just needs to be there
     12224    result = result + recResult;
     12225    for (j = 1; j<=size(recResult); j++)
     12226    {//concatenate i to all elements in recResult
     12227      if (size(recResult[j]) + 1 <= upToDeg)
     12228      {//only concatenate if uptodeg is not exceeded
     12229    recResult[j] = intvec(i, recResult[j]);
     12230      }//only concatenate if uptodeg is not exceeded
     12231    }//concatenate i to all elements in recResult
     12232    result = result + recResult;
     12233  }//For each coordinate not equal to zero, do a recursive call and concatenate
     12234  return(delete_duplicates_noteval_and_sort(result));
     12235}//wordsWithMaxAppearance
     12236
     12237static proc test_wordsWithMaxAppearance()
     12238{//test_wordsWithMaxAppearance
     12239  int result = 1;
     12240  //Test 1: 0 intvec
     12241  list expected = list();
     12242  list obtained = wordsWithMaxAppearance(intvec(0), 4);
     12243  if (!isListEqual(expected,obtained))
     12244  {
     12245    print("Test 1 for wordsWithMaxAppearance failed.");
     12246    print("Expected:\n");
     12247    print(expected);
     12248    print("obtained:\n");
     12249    print(obtained);
     12250    result = 0;
     12251  }
     12252  //Test 2: 0 multiple entries intvec
     12253  expected = list();
     12254  obtained = wordsWithMaxAppearance(0:10, 4);
     12255  if (!isListEqual(expected,obtained))
     12256  {
     12257    print("Test 2 for wordsWithMaxAppearance failed.");
     12258    print("Expected:\n");
     12259    print(expected);
     12260    print("obtained:\n");
     12261    print(obtained);
     12262    result = 0;
     12263  }
     12264  //Test 3: Single letter
     12265  expected = list(intvec(4));
     12266  obtained = wordsWithMaxAppearance(intvec(0, 0, 0, 1, 0, 0), 4);
     12267  if (!isListEqual(expected,obtained))
     12268  {
     12269    print("Test 3 for wordsWithMaxAppearance failed.");
     12270    print("Expected:\n");
     12271    print(expected);
     12272    print("obtained:\n");
     12273    print(obtained);
     12274    result = 0;
     12275  }
     12276  //Test 4: Two letters, one appearance
     12277  expected = list(intvec(3), intvec(3,4), intvec(4), intvec(4,3));
     12278  obtained = wordsWithMaxAppearance(intvec(0, 0, 1, 1, 0, 0), 4);
     12279  if (!isListEqual(expected,obtained))
     12280  {
     12281    print("Test 4 for wordsWithMaxAppearance failed.");
     12282    print("Expected:\n");
     12283    print(expected);
     12284    print("obtained:\n");
     12285    print(obtained);
     12286    result = 0;
     12287  }
     12288  //Test 5: One letter, degree 2
     12289  expected = list(intvec(4), intvec(4,4));
     12290  obtained = wordsWithMaxAppearance(intvec(0, 0, 0, 2, 0, 0), 4);
     12291  if (!isListEqual(expected,obtained))
     12292  {
     12293    print("Test 5 for wordsWithMaxAppearance failed.");
     12294    print("Expected:\n");
     12295    print(expected);
     12296    print("obtained:\n");
     12297    print(obtained);
     12298    result = 0;
     12299  }
     12300  //Test 6: random example
     12301  expected = list(intvec(1),
     12302          intvec(1,4),
     12303          intvec(1,4,4),
     12304          intvec(4),
     12305          intvec(4,1),
     12306          intvec(4,1,4),
     12307          intvec(4,4),
     12308          intvec(4,4,1));
     12309  obtained = wordsWithMaxAppearance(intvec(1, 0, 0, 2, 0, 0), 4);
     12310  if (!isListEqual(expected,obtained))
     12311  {
     12312    print("Test 6 for wordsWithMaxAppearance failed.");
     12313    print("Expected:\n");
     12314    print(expected);
     12315    print("obtained:\n");
     12316    print(obtained);
     12317    result = 0;
     12318  }
     12319  return(result);
     12320}//test_wordsWithMaxAppearance
     12321
     12322static proc monsSmallerThanLetterPlace(poly e, intvec maxDegInCoordinates, int upToDeg)
     12323"USAGE: monsSmallerThanLetterPlace(e, maxDegInCoordinates); e is a
     12324monomial.
     12325maxDegreeInCoordinates encodes the maximal degree we wish to encounter
     12326in each variable.
     12327RETURN: list
     12328PURPOSE: Computes all monomials in the basering which are degree-wise
     12329smaller than the leading monomial of e.
     12330"{//monsSmallerThanLetterPlace
     12331  int p=printlevel-voice+2;//for dbprint
     12332  int i;
     12333  string dbprintWhitespace = "";
     12334  for (i = 1; i<=voice;i++)
     12335  {dbprintWhitespace = dbprintWhitespace + " ";}
     12336  list result;
     12337  int maxLen = min(sum(maxDegInCoordinates), lpDegBound(basering));
     12338  dbprint(p,dbprintWhitespace + "maxLength: " + string(maxLen));
     12339  list allPossibilities = wordsWithMaxAppearance(maxDegInCoordinates, upToDeg);
     12340  dbprint(p,dbprintWhitespace + "words with max appearance: " + string(allPossibilities));
     12341  poly candidate;
     12342  intvec candidateMaxDeg;
     12343  for (i = 1; i<=size(allPossibilities); i++)
     12344  {//checking for monomials with smaller degree than e
     12345    candidate = iv2lp(allPossibilities[i]);
     12346    dbprint(p,dbprintWhitespace + "Checking candidate: " + string(candidate));
     12347    candidateMaxDeg = getMaxDegreeVecLetterPlace(candidate);
     12348    if (candidate < e && candidateMaxDeg <= maxDegInCoordinates)
     12349    {//successfully found a candidate
     12350      result = insert(result, candidate);
     12351    }//successfully found a candidate
     12352  }//checking for monomials with smaller degree than e
     12353  return(result);
     12354}//monsSmallerThanLetterPlace
     12355
     12356static proc test_monsSmallerThanLetterPlace()
     12357{//test_monsSmallerThanLetterPlace
     12358  int result = 1;
     12359  ring r = 0,(x,y,z,d,w),dp;
     12360  int upToDegBound = 10; // degree bound
     12361  def R = makeLetterplaceRing(upToDegBound);
     12362  setring R;
     12363  //Test 1: 0, 0
     12364  poly input1 = 0;
     12365  intvec input2 = intvec(0,0,0,0,0);
     12366  list expected = list();
     12367  list obtained = monsSmallerThanLetterPlace(input1, input2, upToDegBound);
     12368  if (!isListEqual(expected,obtained))
     12369  {
     12370    print("Test 1 for monsSmallerThanLetterPlace failed.");
     12371    print("Expected:\n");
     12372    print(expected);
     12373    print("obtained:\n");
     12374    print(obtained);
     12375    result = 0;
     12376  }
     12377  //Test 2: 0, nontrivial maxdegree
     12378  input1 = 0;
     12379  input2 = intvec(0,0,2,1,0);
     12380  expected = list();
     12381  obtained = monsSmallerThanLetterPlace(input1, input2, upToDegBound );
     12382  if (!isListEqual(expected,obtained))
     12383  {
     12384    print("Test 2 for monsSmallerThanLetterPlace failed.");
     12385    print("Expected:\n");
     12386    print(expected);
     12387    print("obtained:\n");
     12388    print(obtained);
     12389    result = 0;
     12390  }
     12391  //Test 3: constant, nontrivial maxdegree
     12392  input1 = 7;
     12393  input2 = intvec(0,0,2,1,0);
     12394  expected = list();
     12395  obtained = monsSmallerThanLetterPlace(input1, input2, upToDegBound);
     12396  if (!isListEqual(expected,obtained))
     12397  {
     12398    print("Test 3 for monsSmallerThanLetterPlace failed.");
     12399    print("Expected:\n");
     12400    print(expected);
     12401    print("obtained:\n");
     12402    print(obtained);
     12403    result = 0;
     12404  }
     12405  //Test 4: single variable smallest, nontrivial maxdegree
     12406  input1 = w;
     12407  input2 = intvec(0,0,2,1,0);
     12408  expected = list();
     12409  obtained = monsSmallerThanLetterPlace(input1, input2, upToDegBound );
     12410  if (!isListEqual(expected,obtained))
     12411  {
     12412    print("Test 4 for monsSmallerThanLetterPlace failed.");
     12413    print("Expected:\n");
     12414    print(expected);
     12415    print("obtained:\n");
     12416    print(obtained);
     12417    result = 0;
     12418  }
     12419  //Test 5: single variable highest, nontrivial maxdegree
     12420  input1 = x;
     12421  input2 = intvec(0,0,2,1,0);
     12422  expected = list(d, z);
     12423  obtained = monsSmallerThanLetterPlace(input1, input2, upToDegBound);
     12424  if (!isListEqual(expected,obtained))
     12425  {
     12426    print("Test 5 for monsSmallerThanLetterPlace failed.");
     12427    print("Expected:\n");
     12428    print(expected);
     12429    print("obtained:\n");
     12430    print(obtained);
     12431    result = 0;
     12432  }
     12433  //Test 6: single variable highest, nontrivial maxdegree
     12434  input1 = x;
     12435  input2 = intvec(0,0,2,1,0);
     12436  expected = list(d, z);
     12437  obtained = monsSmallerThanLetterPlace(input1, input2, upToDegBound );
     12438  if (!isListEqual(expected,obtained))
     12439  {
     12440    print("Test 6 for monsSmallerThanLetterPlace failed.");
     12441    print("Expected:\n");
     12442    print(expected);
     12443    print("obtained:\n");
     12444    print(obtained);
     12445    result = 0;
     12446  }
     12447  //Test 7: nontrivial monomial, nontrivial maxdegree
     12448  input1 = z*z*d*d;
     12449  input2 = intvec(0,0,2,1,0);
     12450  expected = list(d*z*z,
     12451          d*z,
     12452          d,
     12453          z*d*z,
     12454          z*d,
     12455          z*z*d,
     12456          z*z,
     12457          z);
     12458  obtained = monsSmallerThanLetterPlace(input1, input2, upToDegBound);
     12459  if (!isListEqual(expected,obtained))
     12460  {
     12461    print("Test 7 for monsSmallerThanLetterPlace failed.");
     12462    print("Expected:\n");
     12463    print(expected);
     12464    print("obtained:\n");
     12465    print(obtained);
     12466    result = 0;
     12467  }
     12468  return(result);
     12469}//test_monsSmallerThanLetterPlace
     12470
     12471static proc ncfactor_letterplace_poly_s(poly h)
     12472"USAGE: ncfactor_letterplace_poly_s(h); h is a polynomial in a
     12473        letterplace algebra
     12474RETURN: list(list)
     12475PURPOSE: Compute all factorizations of h
     12476THEORY: Implements an ansatz-driven factorization method as outlined
     12477by Bell, Heinle and Levandovskyy in \"On Noncommutative Finite
     12478Factorization Domains\".
     12479ASSUME:
     12480- basering is a Letterplace ring
     12481- basefield is such, that factorize() can factor any univariate and
     12482  multivariate commutative polynomial over it.
     12483- h is not a constant or 0.
     12484NOTE:
     12485- Every entry of the output list is a list with factors for one possible factorization.
     12486  The first factor is always a constant (1, if no nontrivial constant
     12487  could be excluded).
     12488SEE ALSO: ncfactor
     12489"{//proc ncfactor_letterplace_poly_s
     12490  int p=printlevel-voice+2;//for dbprint
     12491  int i; int j; int k; int l;
     12492  string dbprintWhitespace = "";
     12493  for (i = 1; i<=voice;i++)
     12494  {dbprintWhitespace = dbprintWhitespace + " ";}
     12495  list result;
     12496  // Producing the set M of the paper (permutation of the variables in
     12497  // the leading monomial)
     12498  dbprint(p,dbprintWhitespace + "Generating the set M.");
     12499  def r = basering;
     12500  intvec maxDegrees = getMaxDegreeVecLetterPlace(h);
     12501  dbprint(p,dbprintWhitespace + "The maximum degrees for each var are given by: "+
     12502      string(maxDegrees));
     12503  intvec expVec = lp2iv(leadmonom(h));
     12504  dbprint(p, dbprintWhitespace + "The leading monomial of h is: " + string(expVec));
     12505  list M;
     12506  for (i = 1; i<size(expVec); i++)
     12507  {//The set M consists of all posssible two-word-combinations of the leading monomial
     12508    M = M + list(list(intvec(expVec[1..i]), intvec(expVec[i+1..size(expVec)])));
     12509  }
     12510  dbprint(p,dbprintWhitespace+"The set M is:");
     12511  dbprint(p, M);
     12512  list list_UpToA;
     12513  list list_UpToB;
     12514  list tempSringlist = ringlist(r);
     12515  list clearSlateRingList = ringlist(r);
     12516  clearSlateRingList[2] = list();
     12517  clearSlateRingList[3] = list();
     12518  if (size(clearSlateRingList) == 6)
     12519  {//some non-commutative relations defined
     12520    clearSlateRingList = delete(clearSlateRingList, 6);
     12521    clearSlateRingList = delete(clearSlateRingList, 5);
     12522  }
     12523  else
     12524  {//It is still possible that we have a 5th entry
     12525    if (size(clearSlateRingList) == 5)
     12526    {
     12527      clearSlateRingList = delete(clearSlateRingList, 5);
     12528    }
     12529  }
     12530  intvec tempIntVec;
     12531  def S; def W; def W2;
     12532  list MEntry;
     12533  poly tempA;
     12534  poly tempB;
     12535  int filterFlag;
     12536  int jj = 1; int sizeOfSol = 1;
     12537  for (i = 1; i<=size(M); i++)
     12538  {//Iterating over all two-combinations in the set M
     12539    MEntry = M[i];
     12540    tempA = iv2lp(MEntry[1]);
     12541    tempB = iv2lp(MEntry[2]);
     12542    dbprint(p,dbprintWhitespace + "The potential leading combination of a is " + string(tempA));
     12543    dbprint(p,dbprintWhitespace + "The potential leading combination of b is " + string(tempB));
     12544    dbprint(p,dbprintWhitespace + "Calculating the monomials smaller than A");
     12545    list_UpToA = monsSmallerThanLetterPlace(leadmonom(tempA),
     12546                        maxDegrees - getMaxDegreeVecLetterPlace(tempB),
     12547                        lpDegBound(r));
     12548    dbprint(p,dbprintWhitespace + "Done, they are: " + string(list_UpToA));
     12549    dbprint(p,dbprintWhitespace + "Calculating the monomials smaller than B");
     12550    list_UpToB = monsSmallerThanLetterPlace(leadmonom(tempB),
     12551                        maxDegrees - getMaxDegreeVecLetterPlace(tempA),
     12552                        lpDegBound(r));
     12553    dbprint(p,dbprintWhitespace + "Done, they are: " + string(list_UpToB));
     12554    //we need to filter out lm(tempA) and lm(tempB);
     12555    for (k = size(list_UpToA); k>0; k--)
     12556    {//removing the leading monomial from list_UpToA and invalid parts
     12557      if (list_UpToA[k] == leadmonom(tempA))
     12558      {//found it
     12559    list_UpToA = delete(list_UpToA,k);
     12560    continue;
     12561      }//found it
     12562      if(sum(getMaxDegreeVecLetterPlace(h)) < sum(getMaxDegreeVecLetterPlace(list_UpToA[k])))
     12563      {//total degree exceeds
     12564    list_UpToA = delete(list_UpToA,k);
     12565    continue;
     12566      }//total degree exceeds
     12567    }//removing the leading monomial from list_UpToA and invalid parts
     12568    for (k = size(list_UpToB); k>0; k--)
     12569    {//removing the leading monomial from list_UpToA
     12570      if (list_UpToB[k] == leadmonom(tempB))
     12571      {//found it
     12572    list_UpToB = delete(list_UpToB,k);
     12573    continue;
     12574      }//found it
     12575      if(sum(getMaxDegreeVecLetterPlace(h)) < sum(getMaxDegreeVecLetterPlace(list_UpToB[k])))
     12576      {//total degree exceeds
     12577    list_UpToB = delete(list_UpToB,k);
     12578    continue;
     12579      }//total degree exceeds
     12580    }//removing the leading monomial from list_UpToA
     12581    if (typeof(tempSringlist[1]) == "int")
     12582    {//No extra parameters in the original ring
     12583      tempSringlist[1] = list(tempSringlist[1],list(), list(list()));
     12584    }//No extra parameters in the original ring
     12585    for (k = size(list_UpToB)+1; k>=0; k--)
     12586    {//fill in b-coefficients as variables
     12587      tempSringlist[1][2] = insert(tempSringlist[1][2],"@b("+string(k)+")");
     12588      clearSlateRingList[2] = insert(clearSlateRingList[2],"@b("+string(k)+")");
     12589    }//fill in b-coefficients as variables
     12590    for (k = size(list_UpToA); k>=0; k--)
     12591    {//fill in a-coefficients as variables
     12592      tempSringlist[1][2] = insert(tempSringlist[1][2],"@a("+string(k)+")");
     12593      clearSlateRingList[2] = insert(clearSlateRingList[2], "@a("+string(k)+")");
     12594    }//fill in a-coefficients as variables
     12595    tempIntVec = 0;
     12596    for(k=1; k<=size(list_UpToA) + size(list_UpToB)+3; k++)
     12597    {tempIntVec[k] = 1;}
     12598    clearSlateRingList[3] =
     12599      list(list("lp",tempIntVec),list("C",intvec(0)));
     12600    if (size(tempSringlist[1][3][1]) == 0)
     12601    {//No values in there before
     12602      tempSringlist[1][3][1] =
     12603    list("lp", tempIntVec);
     12604      tempSringlist[1][4] = ideal(0);
     12605    }//No values in there before
     12606    else
     12607    {//there were other parameters
     12608      tempSringlist[1][3][1][2] = tempSringlist[1][3][1][2], tempIntVec;
     12609    }//there were other parameters
     12610    attrib(tempSringlist, "isLetterplaceRing", attrib(r,"isLetterplaceRing"));
     12611    attrib(tempSringlist, "maxExp", 1);
     12612    if (defined(S)) {
     12613      kill S;
     12614    }
     12615    def S = ring(tempSringlist); setring S;
     12616    attrib(S, "uptodeg", lpDegBound(r));
     12617    attrib(S, "isLetterplaceRing", attrib(r,"isLetterplaceRing"));
     12618    dbprint(p,dbprintWhitespace+"Done generating ring S:");
     12619    dbprint(p,dbprintWhitespace+string(S));
     12620    dbprint(p,dbprintWhitespace+"Generate the ansatz.");
     12621    dbprint(p,dbprintWhitespace+"The new ring is: " + string(basering));
     12622    setring(S);
     12623    poly h = imap(r, h);
     12624    poly a = imap(r,tempA);
     12625    poly b = imap(r,tempB);
     12626    if (!defined(list_UpToA))
     12627    {list list_UpToA = imap(r,list_UpToA);}
     12628    if (!defined(list_UpToB))
     12629    {list list_UpToB = imap(r,list_UpToB);}
     12630    b = b*par(size(list_UpToA)+size(list_UpToB)+3);
     12631    for (k = size(list_UpToB); k>0; k--)
     12632    {b = b + par(size(list_UpToA)+2+k)*list_UpToB[k];}
     12633    b = b + par(size(list_UpToA) + 2);
     12634    dbprint(p,dbprintWhitespace + "The ansatz for b is " + string(b));
     12635    for (k = size(list_UpToA); k>0; k--)
     12636    {a = a + par(k+1)*list_UpToA[k];}
     12637    a = a + par(1);
     12638    dbprint(p,dbprintWhitespace + "The ansatz for a is " + string(a));
     12639    poly ansatzPoly = a*b - h;
     12640    dbprint(p,dbprintWhitespace + "The ansatzpoly is " + string(ansatzPoly));
     12641    ideal listOfRelations = coeffs(ansatzPoly, var(1));
     12642    for (k = 2;k<=nvars(r);k++)
     12643    {listOfRelations = coeffs(listOfRelations,var(k));}
     12644    setring(r);
     12645    W = ring(clearSlateRingList); // comm ring of coeffs of the ansatz
     12646    W2 = W+r; // first coeffs, then orig letterplace vars
     12647    setring(W);
     12648    ideal listOfRelations = imap(S, listOfRelations);
     12649    option(redSB);
     12650    option(redTail);
     12651    dbprint(p,dbprintWhitespace + "Calculating the solutions of:");
     12652    dbprint(p,listOfRelations);
     12653    if (!defined(sol))
     12654    {
     12655      def sol = facstd(listOfRelations);
     12656    }
     12657    else
     12658    {
     12659      sol = facstd(listOfRelations);
     12660    }
     12661    dbprint(p,dbprintWhitespace + "Filtering the solutions that are " +
     12662      "not in the base-field.");
     12663    for(k = 1; k <=size(sol);k++)
     12664    {//filtering what is not in the ground-field
     12665      filterFlag = 0;
     12666      for (l=1; l<=size(sol[k]); l++)
     12667      {
     12668        if ((deg(sol[k][l])>1) or ((sol[k][l]/leadcoef(sol[k][l])==1)))
     12669        {//Not in the ground field
     12670          filterFlag = 1;
     12671          break;
     12672        }//Not in the ground field
     12673      }
     12674      if (filterFlag or (vdim(std(sol[k]))<0))
     12675      {//In this case, we can delete that whole entry and move on
     12676        sol= delete(sol,k);
     12677        continue;
     12678      }//In this case, we can delete that whole entry and move on
     12679    }//filtering what is not in the ground-field
     12680    dbprint(p,dbprintWhitespace + "Solutions for the coefficients are:");
     12681    dbprint(p,sol);
     12682    setring S;
     12683    ideal a_coef; jj=1;
     12684    while (a[jj]!=0)
     12685    {
     12686      a_coef[jj]=leadcoef(a[jj]);
     12687      jj++;
     12688    }
     12689    ideal b_coef;
     12690    jj=1;
     12691    while (b[jj]!=0)
     12692    {
     12693      b_coef[jj]=leadcoef(b[jj]);
     12694      jj++;
     12695    }
     12696    setring W; //W = ring of coeffs_variables
     12697    ideal a_coef = imap(S, a_coef);
     12698    ideal b_coef = imap(S, b_coef);
     12699    sizeOfSol = size(sol);
     12700    if (!defined(tempResultCoeffs)) {
     12701      list tempResultCoeffs;
     12702    }
     12703    tempResultCoeffs = list();
     12704    for (k=1; k<=sizeOfSol; k++)
     12705    {
     12706      sol[k] = std(sol[k]);
     12707      tempResultCoeffs = insert(tempResultCoeffs, list(
     12708      NF(a_coef, sol[k]), // element-wise
     12709      NF(b_coef, sol[k]))); // element-wise
     12710    }
     12711    setring S;
     12712    if(!defined(tempResultCoeffs)) {
     12713      list tempResultCoeffs = imap(W, tempResultCoeffs);
     12714    }
     12715    if (!(defined(tempResult)))
     12716    {list tempResult;}
     12717    poly ared;
     12718    poly bred;
     12719    for (k = 1; k <= size(tempResultCoeffs); k++)
     12720    {//Creating the resulting polynomials with the respective coeffs.
     12721      jj=1;
     12722      ared=0;
     12723      bred=0;
     12724      while (a[jj]!=0)
     12725      {//Going through all monomials of a
     12726        ared= ared+leadcoef(tempResultCoeffs[k][1][jj])*leadmonom(a[jj]);
     12727        jj++;
     12728      }//Going through all monomials of a
     12729      jj=1;
     12730      while (b[jj]!=0)
     12731      {//Going through all monomials of b
     12732        bred= bred+leadcoef(tempResultCoeffs[k][2][jj])*leadmonom(b[jj]);
     12733        jj++;
     12734      }//Going through all monomials of b
     12735      tempResult = insert(tempResult,list(ared,bred));
     12736    }//Creating the resulting polynomials with the respective coeffs.
     12737    setring(r);
     12738    if (!defined(tempResult))
     12739    {def tempResult = imap(S,tempResult);}
     12740    for (k=1; k<= size(tempResult);k++)
     12741    {
     12742      result = insert(result,tempResult[k]);
     12743    }
     12744    kill tempResult;
     12745    clearSlateRingList[2] = list();
     12746    clearSlateRingList[3] = list();
     12747    tempSringlist = ringlist(r);
     12748  }//Iterating over all two-combinations in the set M
     12749  if(size(result)==0)
     12750  {
     12751    return (list(list(h)));
     12752  }
     12753  return(delete_duplicates_noteval_and_sort(result));
     12754}//proc ncfactor_letterplace_poly_s
     12755
     12756
     12757static proc test_ncfactor_letterplace_poly_s()
     12758"Tests the basic functionality."
     12759{//proc test_ncfactor_letterplace_poly_s()
     12760  int result = 1;
     12761  ring r = 0,(x,y,z),dp;
     12762  int d =4; // degree bound
     12763  def R = makeLetterplaceRing(d);
     12764  setring R;
     12765  //TEST 1: power of 2
     12766  poly f1 = (x + y)*(x + y);
     12767  list obtained = ncfactor_letterplace_poly_s(f1);
     12768  list expected = sortFactorizations(
     12769    list(
     12770      list(x + y, x + y)
     12771      ));
     12772  if (!isListEqual(expected,obtained))
     12773  {
     12774    print("Test 1 for ncfactor_letterplace_poly_s failed.");
     12775    print("Expected:\n");
     12776    print(expected);
     12777    print("obtained:\n");
     12778    print(obtained);
     12779    result = 0;
     12780  }
     12781  //TEST 2: Monomial
     12782  f1 = x*z;
     12783  obtained = ncfactor_letterplace_poly_s(f1);
     12784  expected = sortFactorizations(
     12785    list(
     12786      list(x, z)
     12787      ));
     12788  if (!isListEqual(expected,obtained))
     12789  {
     12790    print("Test 2 for ncfactor_letterplace_poly_s failed.");
     12791    print("Expected:\n");
     12792    print(expected);
     12793    print("obtained:\n");
     12794    print(obtained);
     12795    result = 0;
     12796  }
     12797  //TEST 3: Generic product with all variables
     12798  f1 = (x*y + z*x)*(x+ z+ y);
     12799  obtained = ncfactor_letterplace_poly_s(f1);
     12800  expected = sortFactorizations(
     12801    list(
     12802      list(
     12803      x*y + z*x,
     12804      x+ z+ y))
     12805      );
     12806  if (!isListEqual(expected,obtained))
     12807  {
     12808    print("Test 3 for ncfactor_letterplace_poly_s failed.");
     12809    print("Expected:\n");
     12810    print(expected);
     12811    print("obtained:\n");
     12812    print(obtained);
     12813    result = 0;
     12814  }
     12815  //TEST 4: More than one factorization
     12816  f1 = x*y*x + x;
     12817  obtained = ncfactor_letterplace_poly_s(f1);
     12818  expected = sortFactorizations(
     12819    list(
     12820      list(x, y*x + 1),
     12821      list(x*y + 1, x))
     12822      );
     12823  if (!isListEqual(expected,obtained))
     12824  {
     12825    print("Test 4 for ncfactor_letterplace_poly_s failed.");
     12826    print("Expected:\n");
     12827    print(expected);
     12828    print("obtained:\n");
     12829    print(obtained);
     12830    result = 0;
     12831  }
     12832  return (result);
     12833}//proc test_ncfactor_letterplace_poly_s()
     12834
     12835static proc ncfactor_without_opt_letterplace(poly h)
     12836"USAGE: ncfactor_without_opt_letterplace(h); h is a polynomial in a
     12837        letterplace ring over a field k.
     12838RETURN: list(list)
     12839PURPOSE: Compute all factorizations of h, without making any
     12840         sanity-checks
     12841THEORY: Implements an ansatz-driven factorization method as outlined
     12842by Bell, Heinle and Levandovskyy in \"On Noncommutative Finite
     12843Factorization Domains\".
     12844ASSUME:
     12845- the basering is a Letterplace ring.
     12846- k is a ring, such that factorize can factor any univariate and
     12847  multivariate commutative polynomial over k.
     12848- There exists at least one variable in the ring.
     12849- h is not a constant or 0.
     12850NOTE:
     12851- Every entry of the output list is a list with factors for one possible factorization.
     12852  The first factor is always a constant (1, if no nontrivial constant
     12853  could be excluded).
     12854SEE ALSO: facWeyl, facSubWeyl, testNCfac, ncfactor
     12855"{//ncfactor_without_opt_letterplace
     12856  int p = printlevel-voice+2;
     12857  int i; int j; int k; int l;
     12858  string dbprintWhitespace = "";
     12859  for (i = 1; i<=voice;i++)
     12860  {dbprintWhitespace = dbprintWhitespace + " ";}
     12861  number commonCoefficient = content(h);
     12862  poly hath = h/commonCoefficient;
     12863  dbprint(p,dbprintWhitespace + "Calculating all possibilities of h as
     12864 a combination of two factors.");
     12865  list result = ncfactor_letterplace_poly_s(hath);
     12866  result = delete_duplicates_noteval_and_sort(result);
     12867  dbprint(p,dbprintWhitespace + "Done. The result is:");
     12868  dbprint(p,result);
     12869  dbprint(p,dbprintWhitespace+"Recursively check factors for
     12870 irreducibility");
     12871  list recursivetemp;
     12872  int changedSomething;
     12873  for(i = 1; i<=size(result);i++)
     12874  {//recursively factorize factors
     12875    if(size(result[i])>1)
     12876    {//Nontrivial factorization
     12877      for (j=1;j<=size(result[i]);j++)
     12878      {//Factorize every factor
     12879        recursivetemp = ncfactor_letterplace_poly_s(result[i][j]);
     12880        recursivetemp = delete_duplicates_noteval(recursivetemp);
     12881        changedSomething = 0;
     12882        for(k=1; k<=size(recursivetemp);k++)
     12883        {//insert factorized factors
     12884          if(size(recursivetemp[k])>1)
     12885          {//nontrivial
     12886            changedSomething = 1;
     12887            result = insert(result,result[i],i);
     12888            for(l = size(recursivetemp[k]);l>=1;l--)
     12889            {
     12890              result[i+1] = insert(result[i+1],recursivetemp[k][l],j);
     12891            }
     12892            result[i+1] = delete(result[i+1],j);
     12893          }//nontrivial
     12894        }//insert factorized factors
     12895        if (changedSomething)
     12896        {break;}
     12897      }//Factorize every factor
     12898      if(changedSomething)
     12899      {
     12900        result = delete(result,i);
     12901        continue;
     12902      }
     12903    }//Nontrivial factorization
     12904  }//recursively factorize factors
     12905  dbprint(p,dbprintWhitespace +" Done");
     12906  dbprint(p,dbprintWhitespace + "Removing duplicates from the list");
     12907  result = delete_duplicates_noteval_and_sort(result);
     12908  for (i = 1; i<=size(result);i++)
     12909  {//Putting the content everywhere
     12910    result[i] = insert(result[i],commonCoefficient);
     12911  }//Putting the content everywhere
     12912  result = normalizeFactors(result);
     12913  result = delete_duplicates_noteval_and_sort(result);
     12914  dbprint(p,dbprintWhitespace +" Done");
     12915  return(result);
     12916}//ncfactor_without_opt_letterplace
     12917
     12918static proc test_ncfactor_without_opt_letterplace()
     12919"Tests the basic functionality."
     12920{//test_ncfactor_without_opt_letterplace
     12921  int result = 1;
     12922  ring r = 0,(x,y,z),dp;
     12923  int d =4; // degree bound
     12924  def R = makeLetterplaceRing(d);
     12925  setring R;
     12926  //TEST 1: power of 2
     12927  poly f1 = 4*(x + y)*(x + y);
     12928  list obtained = ncfactor_without_opt_letterplace(f1);
     12929  list expected = sortFactorizations(
     12930    list(
     12931      list(number(4), x + y, x + y)
     12932      ));
     12933  if (!isListEqual(expected,obtained))
     12934  {
     12935    print("Test 1 for ncfactor_without_opt_letterplace failed.");
     12936    print("Expected:\n");
     12937    print(expected);
     12938    print("obtained:\n");
     12939    print(obtained);
     12940    result = 0;
     12941  }
     12942  //TEST 2: Monomial
     12943  f1 = 11*x*z;
     12944  obtained = ncfactor_without_opt_letterplace(f1);
     12945  expected = sortFactorizations(
     12946    list(
     12947      list(number(11), x, z)
     12948      ));
     12949  if (!isListEqual(expected,obtained))
     12950  {
     12951    print("Test 2 for ncfactor_without_opt_letterplace failed.");
     12952    print("Expected:\n");
     12953    print(expected);
     12954    print("obtained:\n");
     12955    print(obtained);
     12956    result = 0;
     12957  }
     12958  //TEST 3: Generic product with all variables
     12959  f1 = (2*x*y + 4*z*x)*(x+z+y);
     12960  obtained = ncfactor_without_opt_letterplace(f1);
     12961  expected = sortFactorizations(
     12962    list(
     12963      list(number(2), x*y + 2*z*x, x+z+y)
     12964      ));
     12965  if (!isListEqual(expected,obtained))
     12966  {
     12967    print("Test 3 for ncfactor_without_opt_letterplace failed.");
     12968    print("Expected:\n");
     12969    print(expected);
     12970    print("obtained:\n");
     12971    print(obtained);
     12972    result = 0;
     12973  }
     12974  //TEST 4: More than one factorization
     12975  f1 = 6*x*y*x + 9*x;
     12976  obtained = ncfactor_without_opt_letterplace(f1);
     12977  expected = sortFactorizations(
     12978    list(
     12979      list(number(3), x, 2*y*x + 3),
     12980      list(number(3), 2*x*y + 3, x)
     12981      ));
     12982  if (!isListEqual(expected,obtained))
     12983  {
     12984    print("Test 4 for ncfactor_without_opt_letterplace failed.");
     12985    print("Expected:\n");
     12986    print(expected);
     12987    print("obtained:\n");
     12988    print(obtained);
     12989    result = 0;
     12990  }
     12991  //Test 5: At least one recursive call
     12992  f1 = f1*(x + y);
     12993  obtained = ncfactor_without_opt_letterplace(f1);
     12994  expected = sortFactorizations(
     12995    list(
     12996      list(number(3), x, 2*y*x + 3, x + y),
     12997      list(number(3), 2*x*y + 3, x, x + y)
     12998      ));
     12999  if (!isListEqual(expected,obtained))
     13000  {
     13001    print("Test 5 for ncfactor_without_opt_letterplace failed.");
     13002    print("Expected:\n");
     13003    print(expected);
     13004    print("obtained:\n");
     13005    print(obtained);
     13006    result = 0;
     13007  }
     13008  return (result);
     13009}//test_ncfactor_without_opt_letterplace
     13010
     13011
    1198613012//==================================================
    1198713013// EASY EXAMPLES FOR WEYL ALGEBRA
     
    1221713243
    1221813244 /* recent bufgixes (2) Jan 2018; exs from Viktor and Johannes Hoffmann
    12219  // Original bug: found only (x+d)^2*x as faktorization
     13245 // Original bug: found only (x+d)^2*x as factorization
    1222013246LIB "nctools.lib";
    1222113247ring r = 0,(x,d),dp;
     
    1222513251
    1222613252//======================================================================
    12227 // Hard examples not yet calculatable in feasible amount of time
     13253// Hard examples not yet computable in feasible amount of time
    1222813254//======================================================================
    1222913255
Note: See TracChangeset for help on using the changeset viewer.