Changeset a1f90e4 in git


Ignore:
Timestamp:
Oct 5, 2010, 10:57:28 PM (14 years ago)
Author:
Viktor Levandovskyy <levandov@…>
Branches:
(u'fieker-DuVal', '117eb8c30fc9e991c4decca4832b1d19036c4c65')(u'spielwiese', '8a337797cc4177aa8747d661d5c4214ea2095dac')
Children:
8fc505d5c708f5f706bc2fc4c71e91280d403c61
Parents:
a212fe5c82b2aa44b0cbcac937f3c282fd6d4f10
Message:
*levandov: release preparation factorization in the 1st weyl algebra

git-svn-id: file:///usr/local/Singular/svn/trunk@13397 2c84dea3-7e68-4137-9b89-c4e89433aadc
File:
1 edited

Legend:

Unmodified
Added
Removed
  • Singular/LIB/ncfactor.lib

    ra212fe ra1f90e4  
    1 ////////////////////////////////////////////////////////////
     1///////////////////////////////////////////////////////////
    22version = "$Id$";
    33category="Noncommutative";
     
    77
    88PROCEDURES:
    9 facfirstwa(h);     factorization in the first Weyl algebra
     9facFirstWeyl(h);    factorization in the first Weyl algebra
     10testNCfac(l[,h]);   tests factorizations from a given list for correctness
     11facSubWeyl(h,X,D);  factorization in the first Weyl algebra as a subalgebra
    1012";
    1113
     
    1315LIB "nctools.lib";
    1416LIB "involut.lib";
    15 
    16 //static homogfac(h);         computes one factorization of a homogeneous polynomial h
    17 //static homogfac_all(h);     computes all factorizations of a homogeneous polynomial h
    18 //static facwa(h);            computes all factorizations of an inhomogeneous polynomial h
     17LIB "freegb.lib"; // for isVar
     18
    1919/* ring R = 0,(x,y),Ws(-1,1); */
    2020/* def r = nc_algebra(1,1); */
    2121/* setring(r); */
    2222
    23 ////////////////////////////////////////////////////////////////////////////////////////////////////
     23/////////////////////////////////////////////////////
    2424//==================================================*
    2525//deletes double-entries in a list of factorization
     
    3333     for (i = 1; i<= size(l); i++)
    3434     {//Iterate over the different factorizations
    35           for (j = i+1; j<= size(l); j++)
    36           {//Compare the i'th factorization to the j'th
    37                if (size(l[i])!= size(l[j]))
    38                {//different sizes => not equal
    39                     j++;
    40                     continue;
    41                }//different sizes => not equal
    42                is_equal = 1;
    43                for (k = 1; k <= size(l[i]);k++)
    44                {//Compare every entry
    45                     if (l[i][k]!=l[j][k])
    46                     {
    47                         is_equal = 0;
    48                         break;
    49                     }
    50                }//Compare every entry
    51                if (is_equal == 1)
    52                {//Delete this entry, because there is another equal one int the list
    53                     result = delete(result, i-deleted);
    54                     deleted = deleted+1;
    55                     break;
    56                }//Delete this entry, because there is another equal one int the list
    57           }//Compare the i'th factorization to the j'th
     35          for (j = i+1; j<= size(l); j++)
     36          {//Compare the i'th factorization to the j'th
     37               if (size(l[i])!= size(l[j]))
     38               {//different sizes => not equal
     39                    j++;
     40                    continue;
     41           }//different sizes => not equal
     42           is_equal = 1;
     43           for (k = 1; k <= size(l[i]);k++)
     44           {//Compare every entry
     45            if (l[i][k]!=l[j][k])
     46            {
     47            is_equal = 0;
     48            break;
     49            }
     50           }//Compare every entry
     51           if (is_equal == 1)
     52           {//Delete this entry, because there is another equal one int the list
     53            result = delete(result, i-deleted);
     54            deleted = deleted+1;
     55            break;
     56           }//Delete this entry, because there is another equal one int the list
     57      }//Compare the i'th factorization to the j'th
    5858     }//Iterate over the different factorizations
    5959     return(result);
     
    7171     for (i = 1; i<= size(result); i++)
    7272     {//Iterating over all elements in result
    73           for (j = i+1; j<= size(result); j++)
    74           {//comparing with the other elements
    75                if (product(result[i]) == product(result[j]))
    76                {//There are two equal results; throw away that one with the smaller size
    77                     if (size(result[i])>=size(result[j]))
    78                     {//result[i] has more entries
    79                         result = delete(result,j);
    80                         continue;
    81                     }//result[i] has more entries
    82                     else
    83                     {//result[j] has more entries
    84                         result = delete(result,i);
    85                         i--;
    86                         break;
    87                     }//result[j] has more entries
    88                }//There are two equal results; throw away that one with the smaller size
    89           }//comparing with the other elements
     73      for (j = i+1; j<= size(result); j++)
     74      {//comparing with the other elements
     75           if (product(result[i]) == product(result[j]))
     76           {//There are two equal results; throw away that one with the smaller size
     77            if (size(result[i])>=size(result[j]))
     78            {//result[i] has more entries
     79            result = delete(result,j);
     80            continue;
     81            }//result[i] has more entries
     82            else
     83            {//result[j] has more entries
     84            result = delete(result,i);
     85            i--;
     86            break;
     87            }//result[j] has more entries
     88           }//There are two equal results; throw away that one with the smaller size
     89      }//comparing with the other elements
    9090     }//Iterating over all elements in result
    9191     return(result);
     
    106106     if (nofgl == 0)
    107107     {//g was the empty list
    108           return(result);
     108      return(result);
    109109     }//g was the empty list
    110110     if (nof <= 0)
    111111     {//The user wants to recieve a negative number or no element as a result
    112           return(result);
     112      return(result);
    113113     }//The user wants to recieve a negative number or no element as a result
    114114     if (nofgl == nof)
    115115     {//There are no factors to combine
    116           if (limitcheck(g,limits))
    117           {
    118                result = result + list(g);
    119           }
    120           return(result);
     116      if (limitcheck(g,limits))
     117      {
     118           result = result + list(g);
     119      }
     120      return(result);
    121121     }//There are no factors to combine
    122122     if (nof == 1)
    123123     {//User wants to get just one factor
    124           if (limitcheck(list(product(g)),limits))
    125           {
    126                result = result + list(list(product(g)));
    127           }
    128           return(result);
     124      if (limitcheck(list(product(g)),limits))
     125      {
     126           result = result + list(list(product(g)));
     127      }
     128      return(result);
    129129     }//User wants to get just one factor
    130130     for (i = nof; i > 1; i--)
    131131     {//computing the possibilities that have at least one original factor from g
    132           for (j = i; j>=1; j--)
    133           {//shifting the window of combinable factors to the left
    134                fc = combinekfinlf(list(g[(j)..(j+nofgl - i)]),nof - i + 1,limits); //fc stands for "factors combined"
    135                for (k = 1; k<=size(fc); k++)
    136                {//iterating over the different solutions of the smaller problem
    137                     if (j>1)
    138                     {//There are g_i before the combination
    139                         if (j+nofgl -i < nofgl)
    140                         {//There are g_i after the combination
    141                               temp = list(g[1..(j-1)]) + fc[k] + list(g[(j+nofgl-i+1)..nofgl]);
    142                         }//There are g_i after the combination
    143                         else
    144                         {//There are no g_i after the combination
    145                               temp = list(g[1..(j-1)]) + fc[k];
    146                         }//There are no g_i after the combination
    147                     }//There are g_i before the combination
    148                     if (j==1)
    149                     {//There are no g_i before the combination
    150                         if (j+ nofgl -i <nofgl)
    151                         {//There are g_i after the combination
    152                               temp = fc[k]+ list(g[(j + nofgl - i +1)..nofgl]);
    153                         }//There are g_i after the combination
    154                     }//There are no g_i before the combination
    155                     if (limitcheck(temp,limits))
    156                     {
    157                         result = result + list(temp);
    158                     }
    159                }//iterating over the different solutions of the smaller problem
    160           }//shifting the window of combinable factors to the left
     132      for (j = i; j>=1; j--)
     133      {//shifting the window of combinable factors to the left
     134           fc = combinekfinlf(list(g[(j)..(j+nofgl - i)]),nof - i + 1,limits); //fc stands for "factors combined"
     135           for (k = 1; k<=size(fc); k++)
     136           {//iterating over the different solutions of the smaller problem
     137            if (j>1)
     138            {//There are g_i before the combination
     139            if (j+nofgl -i < nofgl)
     140            {//There are g_i after the combination
     141                  temp = list(g[1..(j-1)]) + fc[k] + list(g[(j+nofgl-i+1)..nofgl]);
     142            }//There are g_i after the combination
     143            else
     144            {//There are no g_i after the combination
     145                  temp = list(g[1..(j-1)]) + fc[k];
     146            }//There are no g_i after the combination
     147            }//There are g_i before the combination
     148            if (j==1)
     149            {//There are no g_i before the combination
     150            if (j+ nofgl -i <nofgl)
     151            {//There are g_i after the combination
     152                  temp = fc[k]+ list(g[(j + nofgl - i +1)..nofgl]);
     153            }//There are g_i after the combination
     154            }//There are no g_i before the combination
     155            if (limitcheck(temp,limits))
     156            {
     157            result = result + list(temp);
     158            }
     159           }//iterating over the different solutions of the smaller problem
     160      }//shifting the window of combinable factors to the left
    161161     }//computing the possibilities that have at least one original factor from g
    162162     for (i = 2; i<=nofgl/nof;i++)
    163163     {//getting the other possible results
    164           result = result + combinekfinlf(list(product(list(g[1..i])))+list(g[(i+1)..nofgl]),nof,limits);
     164      result = result + combinekfinlf(list(product(list(g[1..i])))+list(g[(i+1)..nofgl]),nof,limits);
    165165     }//getting the other possible results
    166166     result = delete_dublicates_noteval(result);
     
    179179     if (size(l1)==0)
    180180     {
    181           return(list());
     181      return(list());
    182182     }
    183183     if (size(l2)==0)
    184184     {
    185           return(list());
     185      return(list());
    186186     }
    187187     if (size(l2)<=size(l1))
    188188     {//l1 will be our g, l2 our f
    189           g = l1;
    190           f = l2;
     189      g = l1;
     190      f = l2;
    191191     }//l1 will be our g, l2 our f
    192192     else
    193193     {//l1 will be our f, l2 our g
    194           g = l2;
    195           f = l1;
     194      g = l2;
     195      f = l1;
    196196     }//l1 will be our f, l2 our g
    197197     def result = combinekfinlf(g,size(f),limits);
    198198     for (i = 1 ; i<= size(result); i++)
    199199     {//Adding the factors of f to every possibility listed in temp
    200           for (j = 1; j<= size(f); j++)
    201           {
    202                result[i][j] = result[i][j]+f[j];
    203           }
    204           if(!limitcheck(result[i],limits))
    205           {
    206                result = delete(result,i);
    207                i--;
    208           }
     200      for (j = 1; j<= size(f); j++)
     201      {
     202           result[i][j] = result[i][j]+f[j];
     203      }
     204      if(!limitcheck(result[i],limits))
     205      {
     206           result = delete(result,i);
     207           i--;
     208      }
    209209     }//Adding the factors of f to every possibility listed in temp
    210210     return(result);
     
    226226     if (size(l1)==0)
    227227     {//the first list is empty
    228           return(list());
     228      return(list());
    229229     }//the first list is empty
    230230     if(size(l2)==0)
    231231     {//the second list is empty
    232           return(list());
     232      return(list());
    233233     }//the second list is empty
    234234     if (size(l2)<=size(l1))
    235235     {//l1 will be our g, l2 our f
    236           g = l1;
    237           f = l2;
     236      g = l1;
     237      f = l2;
    238238     }//l1 will be our g, l2 our f
    239239     else
    240240     {//l1 will be our f, l2 our g
    241           g = l2;
    242           f = l1;
     241      g = l2;
     242      f = l1;
    243243     }//l1 will be our f, l2 our g
    244244     list M;
    245245     for (i = 2; i<size(f); i++)
    246246     {//finding common factors of f and g...
    247         for (j=2; j<size(g);j++)
    248         {//... with g
    249              if (f[i] == g[j])
    250              {//we have an equal pair
    251                   M = M + list(list(i,j));
    252              }//we have an equal pair
    253         }//... with g
     247    for (j=2; j<size(g);j++)
     248    {//... with g
     249         if (f[i] == g[j])
     250         {//we have an equal pair
     251          M = M + list(list(i,j));
     252         }//we have an equal pair
     253    }//... with g
    254254     }//finding common factors of f and g...
    255255     if (g[1]==f[1])
    256256     {//Checking for the first elements to be equal
    257           M = M + list(list(1,1));
     257      M = M + list(list(1,1));
    258258     }//Checking for the first elements to be equal
    259259     if (g[size(g)]==f[size(f)])
    260260     {//Checking for the last elements to be equal
    261           M = M + list(list(size(f),size(g)));
     261      M = M + list(list(size(f),size(g)));
    262262     }//Checking for the last elements to be equal
    263263     list result;//= list(list());
    264264     while(size(M)>0)
    265265     {//set of equal pairs is not empty
    266           temp = M[1];
    267           temppos = 1;
    268           for (i = 2; i<=size(M); i++)
    269           {//finding the minimal element of M
    270                if (M[i][1]<=temp[1])
    271                {//a possible candidate that is smaller than temp could have been found
    272                     if (M[i][1]==temp[1])
    273                     {//In this case we must look at the second number
    274                         if (M[i][2]< temp[2])
    275                         {//the candidate is smaller
    276                               temp = M[i];
    277                               temppos = i;
    278                         }//the candidate is smaller
    279                     }//In this case we must look at the second number
    280                     else
    281                     {//The candidate is definately smaller
    282                         temp = M[i];
    283                         temppos = i;
    284                     }//The candidate is definately smaller
    285                }//a possible candidate that is smaller than temp could have been found
    286           }//finding the minimal element of M
    287           M = delete(M, temppos);
    288           if(temp[1]>1)
    289           {//There are factors to combine before the equal factor
    290                if (temp[1]<size(f))
    291                {//The most common case
    292                     //first the combinations ignoring common factors
    293                     pre = merge_icf(list(f[1..(temp[1]-1)]),list(g[1..(temp[2]-1)]),limits);
    294                     post = merge_icf(list(f[(temp[1]+1)..size(f)]),list(g[(temp[2]+1..size(g))]),limits);
    295                     for (i = 1; i <= size(pre); i++)
    296                     {//all possible pre's...
    297                         for (j = 1; j<= size(post); j++)
    298                         {//...combined with all possible post's
    299                               candidate = pre[i]+list(f[temp[1]])+post[j];
    300                               if (limitcheck(candidate,limits))
    301                               {
    302                                    result = result + list(candidate);
    303                               }
    304                         }//...combined with all possible post's
    305                     }//all possible pre's...
    306                     //Now the combinations with respect to common factors
    307                     post = merge_cf(list(f[(temp[1]+1)..size(f)]),list(g[(temp[2]+1..size(g))]),limits);
    308                     if (size(post)>0)
    309                     {//There are factors to combine
    310                         for (i = 1; i <= size(pre); i++)
    311                         {//all possible pre's...
    312                               for (j = 1; j<= size(post); j++)
    313                               {//...combined with all possible post's
    314                                    candidate= pre[i]+list(f[temp[1]])+post[j];
    315                                    if (limitcheck(candidate,limits))
    316                                    {
    317                                         result = result + list(candidate);
    318                                    }
    319                               }//...combined with all possible post's
    320                         }//all possible pre's...
    321                     }//There are factors to combine
    322                }//The most common case
    323                else
    324                {//the last factor is the common one
    325                     pre = merge_icf(list(f[1..(temp[1]-1)]),list(g[1..(temp[2]-1)]),limits);
    326                     for (i = 1; i<= size(pre); i++)
    327                     {//iterating over the possible pre-factors
    328                         candidate = pre[i]+list(f[temp[1]]);
    329                         if (limitcheck(candidate,limits))
    330                         {
    331                               result = result + list(candidate);
    332                         }
    333                     }//iterating over the possible pre-factors
    334                }//the last factor is the common one
    335           }//There are factors to combine before the equal factor
    336           else
    337           {//There are no factors to combine before the equal factor
    338                if (temp[1]<size(f))
    339                {//Just a check for security
    340                     //first without common factors
    341                     post=merge_icf(list(f[(temp[1]+1)..size(f)]),list(g[(temp[2]+1..size(g))]),limits);
    342                     for (i = 1; i<=size(post); i++)
    343                     {
    344                         candidate = list(f[temp[1]])+post[i];
    345                         if (limitcheck(candidate,limits))
    346                         {
    347                               result = result + list(candidate);
    348                         }
    349                     }
    350                     //Now with common factors
    351                     post = merge_cf(list(f[(temp[1]+1)..size(f)]),list(g[(temp[2]+1..size(g))]),limits);
    352                     if(size(post)>0)
    353                     {//we could find other combinations
    354                         for (i = 1; i<=size(post); i++)
    355                         {
    356                               candidate = list(f[temp[1]])+post[i];
    357                               if (limitcheck(candidate,limits))
    358                               {
    359                                    result = result + list(candidate);
    360                               }
    361                         }
    362                     }//we could find other combinations
    363                }//Just a check for security
    364           }//There are no factors to combine before the equal factor
     266      temp = M[1];
     267      temppos = 1;
     268      for (i = 2; i<=size(M); i++)
     269      {//finding the minimal element of M
     270           if (M[i][1]<=temp[1])
     271           {//a possible candidate that is smaller than temp could have been found
     272            if (M[i][1]==temp[1])
     273            {//In this case we must look at the second number
     274            if (M[i][2]< temp[2])
     275                    {//the candidate is smaller
     276                              temp = M[i];
     277                              temppos = i;
     278                        }//the candidate is smaller
     279                    }//In this case we must look at the second number
     280                    else
     281                    {//The candidate is definately smaller
     282                        temp = M[i];
     283                        temppos = i;
     284                    }//The candidate is definately smaller
     285               }//a possible candidate that is smaller than temp could have been found
     286          }//finding the minimal element of M
     287          M = delete(M, temppos);
     288          if(temp[1]>1)
     289          {//There are factors to combine before the equal factor
     290               if (temp[1]<size(f))
     291               {//The most common case
     292                    //first the combinations ignoring common factors
     293                    pre = merge_icf(list(f[1..(temp[1]-1)]),list(g[1..(temp[2]-1)]),limits);
     294                    post = merge_icf(list(f[(temp[1]+1)..size(f)]),list(g[(temp[2]+1..size(g))]),limits);
     295                    for (i = 1; i <= size(pre); i++)
     296                    {//all possible pre's...
     297                        for (j = 1; j<= size(post); j++)
     298                        {//...combined with all possible post's
     299                              candidate = pre[i]+list(f[temp[1]])+post[j];
     300                              if (limitcheck(candidate,limits))
     301                              {
     302                                   result = result + list(candidate);
     303                              }
     304                        }//...combined with all possible post's
     305                    }//all possible pre's...
     306                    //Now the combinations with respect to common factors
     307                    post = merge_cf(list(f[(temp[1]+1)..size(f)]),list(g[(temp[2]+1..size(g))]),limits);
     308                    if (size(post)>0)
     309                    {//There are factors to combine
     310                        for (i = 1; i <= size(pre); i++)
     311                        {//all possible pre's...
     312                              for (j = 1; j<= size(post); j++)
     313                              {//...combined with all possible post's
     314                                   candidate= pre[i]+list(f[temp[1]])+post[j];
     315                                   if (limitcheck(candidate,limits))
     316                                   {
     317                                        result = result + list(candidate);
     318                                   }
     319                              }//...combined with all possible post's
     320                        }//all possible pre's...
     321                    }//There are factors to combine
     322               }//The most common case
     323               else
     324               {//the last factor is the common one
     325                    pre = merge_icf(list(f[1..(temp[1]-1)]),list(g[1..(temp[2]-1)]),limits);
     326                    for (i = 1; i<= size(pre); i++)
     327                    {//iterating over the possible pre-factors
     328                        candidate = pre[i]+list(f[temp[1]]);
     329                        if (limitcheck(candidate,limits))
     330                        {
     331                              result = result + list(candidate);
     332                        }
     333                    }//iterating over the possible pre-factors
     334               }//the last factor is the common one
     335          }//There are factors to combine before the equal factor
     336          else
     337          {//There are no factors to combine before the equal factor
     338               if (temp[1]<size(f))
     339               {//Just a check for security
     340                    //first without common factors
     341                    post=merge_icf(list(f[(temp[1]+1)..size(f)]),list(g[(temp[2]+1..size(g))]),limits);
     342                    for (i = 1; i<=size(post); i++)
     343                    {
     344                        candidate = list(f[temp[1]])+post[i];
     345                        if (limitcheck(candidate,limits))
     346                        {
     347                              result = result + list(candidate);
     348                        }
     349                    }
     350                    //Now with common factors
     351                    post = merge_cf(list(f[(temp[1]+1)..size(f)]),list(g[(temp[2]+1..size(g))]),limits);
     352                    if(size(post)>0)
     353                    {//we could find other combinations
     354                        for (i = 1; i<=size(post); i++)
     355                        {
     356                              candidate = list(f[temp[1]])+post[i];
     357                              if (limitcheck(candidate,limits))
     358                              {
     359                                   result = result + list(candidate);
     360                              }
     361                        }
     362                    }//we could find other combinations
     363               }//Just a check for security
     364          }//There are no factors to combine before the equal factor
    365365     }//set of equal pairs is not empty
    366366     return(result);
     
    379379     if (size(l2)<=size(l1))
    380380     {//l1 will be our g, l2 our f
    381           g = l1;
    382           f = l2;
     381          g = l1;
     382          f = l2;
    383383     }//l1 will be our g, l2 our f
    384384     else
    385385     {//l1 will be our f, l2 our g
    386           g = l2;
    387           f = l1;
     386          g = l2;
     387          f = l1;
    388388     }//l1 will be our f, l2 our g
    389389     list result;
    390390     for (l = size(f); l>=1; l--)
    391391     {//all possibilities to combine the factors of f
    392           F = combinekfinlf(f,l,limits);
    393           for (k = 1; k<= size(F);k++)
     392          F = combinekfinlf(f,l,limits);
     393          for (k = 1; k<= size(F);k++)
    394394          {//for all possibilities of combinations of the factors of f
    395                result = result + merge_cf(F[k],g,limits);
    396                result = result + merge_icf(F[k],g,limits);
    397           }//for all possibilities of combinations of the factors of f
     395               result = result + merge_cf(F[k],g,limits);
     396               result = result + merge_icf(F[k],g,limits);
     397          }//for all possibilities of combinations of the factors of f
    398398     }//all possibilities to combine the factors of f
    399399     return(result);
     
    408408     if (size(limits)!=3)
    409409     {//check the input
    410           return(0);
     410          return(0);
    411411     }//check the input
    412412     if(size(g)==0)
    413413     {
    414           return(0);
     414          return(0);
    415415     }
    416416     def prod = product(g);
     
    418418     for (i = 1; i<=size(limg);i++)
    419419     {//the final check
    420           if(limg[i]>limits[i])
    421           {
    422                return(0);
    423           }
     420          if(limg[i]>limits[i])
     421          {
     422               return(0);
     423          }
    424424     }//the final check
    425425     return(1);
     
    431431//in the first Weyl Algebra
    432432static proc homogfac(poly h)
    433 "USAGE: homogfac(h); h is a homogeneous polynomial in the first Weyl algebra with respect to the weight vector [-1,1]
     433"USAGE: homogfac(h); h is a homogeneous polynomial in the
     434        first Weyl algebra with respect to the weight vector [-1,1]
    434435RETURN: list
    435 PURPOSE: Computes a factorization of a homogeneous polynomial h with respect to the weight vector [-1,1] in the first Weyl algebra
    436 THEORY: @code{homogfac} returns a list with a factorization of the given, homogeneous polynomial. If the degree of the polynomial is k with k positive, the last k entries in the output list are the second variable. If k is positive, the last k entries will be x. The other entries will be irreducible polynomials of degree zero or 1 resp. -1.
    437 EXAMPLE: example homogfac; shows examples
     436PURPOSE: Computes a factorization of a homogeneous polynomial h with
     437         respect to the weight vector [-1,1] in the first Weyl algebra
     438THEORY: @code{homogfac} returns a list with a factorization of the given,
     439        [-1,1]-homogeneous polynomial. If the degree of the polynomial is k with
     440        k positive, the last k entries in the output list are the second
     441        variable. If k is positive, the last k entries will be x. The other
     442        entries will be irreducible polynomials of degree zero or 1 resp. -1.
    438443SEE ALSO: homogfac_all
    439444"{//proc homogfac
     
    443448     if (!homogwithorder(h,intvec(-1,1)))
    444449     {//The given polynomial is not homogeneous
    445           return(list());
     450          return(list());
    446451     }//The given polynomial is not homogeneous
    447452     if (h==0)
    448453     {
    449           return(list(0));
     454          return(list(0));
    450455     }
    451456     list result;
     
    453458     if (m!=0)
    454459     {//The degree is not zero
    455           if (m <0)
    456           {//There are more x than y
    457                hath = lift(var(1)^(-m),h)[1,1];
    458                for (i = 1; i<=-m; i++)
    459                {
    460                     result = result + list(var(1));
    461                }
    462           }//There are more x than y
    463           else
    464           {//There are more y than x
    465                hath = lift(var(2)^m,h)[1,1];
    466                for (i = 1; i<=m;i++)
    467                {
    468                     result = result + list(var(2));
    469                }
    470           }//There are more y than x
     460          if (m <0)
     461          {//There are more x than y
     462               hath = lift(var(1)^(-m),h)[1,1];
     463               for (i = 1; i<=-m; i++)
     464               {
     465                    result = result + list(var(1));
     466               }
     467          }//There are more x than y
     468          else
     469          {//There are more y than x
     470               hath = lift(var(2)^m,h)[1,1];
     471               for (i = 1; i<=m;i++)
     472               {
     473                    result = result + list(var(2));
     474               }
     475          }//There are more y than x
    471476     }//The degree is not zero
    472477     else
    473478     {//The degree is zero
    474           hath = h;
     479          hath = h;
    475480     }//The degree is zero
    476481     //beginning to transform x^i*y^i in teta(teta-1)...(teta-i+1)
     
    478483     for(i = 1; i<=size(hath);i++)
    479484     {//Putting the monomials in a list
    480           mons = mons+list(hath[i]);
     485          mons = mons+list(hath[i]);
    481486     }//Putting the monomials in a list
    482487     ring tempRing = 0,(x,y,teta),dp;
     
    487492     for (i = 1; i<=size(mons);i++)
    488493     {//transforming the monomials as monomials in theta
    489           entry = leadcoef(mons[i]);
    490           for (j = 0; j<leadexp(mons[i])[2];j++)
    491           {
    492                entry = entry * (teta-j);
    493           }
    494           mons[i] = entry;
     494          entry = leadcoef(mons[i]);
     495          for (j = 0; j<leadexp(mons[i])[2];j++)
     496          {
     497               entry = entry * (teta-j);
     498          }
     499          mons[i] = entry;
    495500     }//transforming the monomials as monomials in theta
    496501     list azeroresult = factorize(sum(mons));
     
    498503     for (i = 1; i<=size(azeroresult[1]);i++)
    499504     {//rewrite the result of the commutative factorization
    500           for (j = 1; j <= azeroresult[2][i];j++)
    501           {
    502                azeroresult_return_form = azeroresult_return_form + list(azeroresult[1][i]);
    503           }
     505          for (j = 1; j <= azeroresult[2][i];j++)
     506          {
     507               azeroresult_return_form = azeroresult_return_form + list(azeroresult[1][i]);
     508          }
    504509     }//rewrite the result of the commutative factorization
    505510     setring(r);
     
    508513     for (i = 1; i<=size(tempresult);i++)
    509514     {//factorizations of theta resp. theta +1
    510           if(tempresult[i]==var(1)*var(2))
    511           {
    512                tempresult = insert(tempresult,var(1),i-1);
    513                i++;
    514                tempresult[i]=var(2);
    515           }
    516           if(tempresult[i]==var(2)*var(1))
    517           {
    518                tempresult = insert(tempresult,var(2),i-1);
    519                i++;
    520                tempresult[i]=var(1);
    521           }
     515          if(tempresult[i]==var(1)*var(2))
     516          {
     517               tempresult = insert(tempresult,var(1),i-1);
     518               i++;
     519               tempresult[i]=var(2);
     520          }
     521          if(tempresult[i]==var(2)*var(1))
     522          {
     523               tempresult = insert(tempresult,var(2),i-1);
     524               i++;
     525               tempresult[i]=var(1);
     526          }
    522527     }//factorizations of theta resp. theta +1
    523528     result = tempresult+result;
     
    537542//Computes all possible homogeneous factorizations
    538543static proc homogfac_all(poly h)
    539 "USAGE: homogfac_all(h); h is a homogeneous polynomial in the first Weyl algebra with respect to the weight vector [-1,1]
     544"USAGE: homogfac_all(h); h is a homogeneous polynomial in the first Weyl algebra
     545        with respect to the weight vector [-1,1]
    540546RETURN: list
    541 PURPOSE: Computes all factorizations of a homogeneous polynomial h with respect to the weight vector [-1,1] in the first Weyl algebra
    542 THEORY: @code{homogfac} returns a list with all factorization of the given, homogeneous polynomial. It uses the output of homogfac and permutes its entries with respect to the commutation rule. Furthermore, if a factor of degree zero is irreducible in K[\theta], but reducible in the first Weyl algebra, the permutations of this element with the other entries will also be computed.
    543 EXAMPLE: example homogfac; shows examples
     547PURPOSE: Computes all factorizations of a homogeneous polynomial h with respect
     548         to the weight vector [-1,1] in the first Weyl algebra
     549THEORY: @code{homogfac} returns a list with all factorization of the given,
     550        homogeneous polynomial. It uses the output of homogfac and permutes
     551        its entries with respect to the commutation rule. Furthermore, if a
     552        factor of degree zero is irreducible in K[\theta], but reducible in
     553        the first Weyl algebra, the permutations of this element with the other
     554        entries will also be computed.
    544555SEE ALSO: homogfac
    545556"{//proc HomogFacAll
    546557     if (deg(h,intvec(1,1)) <= 0 )
    547558     {//h is a constant
    548           return(list(list(h)));
     559          return(list(list(h)));
    549560     }//h is a constant
    550561     def r = basering;
     
    555566     if (size(one_hom_fac) == 0)
    556567     {//there is no homogeneous factorization or the polynomial was not homogeneous
    557           return(list());
     568          return(list());
    558569     }//there is no homogeneous factorization or the polynomial was not homogeneous
    559570     //divide list in A0-Part and a list of x's resp. y's
     
    566577     if (absValue(deg(h,intvec(-1,1)))<size(one_hom_fac)-1)
    567578     {//There is a nontrivial A_0-part
    568           list_azero = one_hom_fac[2..(size(one_hom_fac)-absValue(deg(h,intvec(-1,1))))];
    569           is_list_azero_empty = 0;
     579          list_azero = one_hom_fac[2..(size(one_hom_fac)-absValue(deg(h,intvec(-1,1))))];
     580          is_list_azero_empty = 0;
    570581     }//There is a nontrivial A_0 part
    571582     for (i = 1; i<=size(list_azero)-1;i++)
    572583     {//in homogfac, we factorized theta, and this will be made undone
    573           if (list_azero[i] == var(1))
    574           {
    575                if (list_azero[i+1]==var(2))
    576                {
    577                     list_azero[i] = var(1)*var(2);
    578                     list_azero = delete(list_azero,i+1);
    579                }
    580           }
    581           if (list_azero[i] == var(2))
    582           {
    583                if (list_azero[i+1]==var(1))
    584                {
    585                     list_azero[i] = var(2)*var(1);
    586                     list_azero = delete(list_azero,i+1);
    587                }
    588           }
     584          if (list_azero[i] == var(1))
     585          {
     586               if (list_azero[i+1]==var(2))
     587               {
     588                    list_azero[i] = var(1)*var(2);
     589                    list_azero = delete(list_azero,i+1);
     590               }
     591          }
     592          if (list_azero[i] == var(2))
     593          {
     594               if (list_azero[i+1]==var(1))
     595               {
     596                    list_azero[i] = var(2)*var(1);
     597                    list_azero = delete(list_azero,i+1);
     598               }
     599          }
    589600     }//in homogfac, we factorized theta, and this will be made undone
    590601     if(deg(h,intvec(-1,1))!=0)
    591602     {//list_not_azero is not empty
    592           list_not_azero = one_hom_fac[(size(one_hom_fac)-absValue(deg(h,intvec(-1,1)))+1)..size(one_hom_fac)];
    593           is_list_not_azero_empty = 0;
     603          list_not_azero = one_hom_fac[(size(one_hom_fac)-absValue(deg(h,intvec(-1,1)))+1)..size(one_hom_fac)];
     604          is_list_not_azero_empty = 0;
    594605     }//list_not_azero is not empty
    595606     //Map list_azero in K[theta]
     
    600611     if(!is_list_not_azero_empty)
    601612     {//Mapping in Singular is only possible, if the list before contained at least one element of the other ring
    602           list list_not_azero = thetamap(list_not_azero);
     613          list list_not_azero = thetamap(list_not_azero);
    603614     }//Mapping in Singular is only possible, if the list before contained at least one element of the other ring
    604615     if(!is_list_azero_empty)
    605616     {//Mapping in Singular is only possible, if the list before contained at least one element of the other ring
    606           list list_azero= thetamap(list_azero);
     617          list list_azero= thetamap(list_azero);
    607618     }//Mapping in Singular is only possible, if the list before contained at least one element of the other ring
    608619     list k_factor = thetamap(k_factor);
     
    610621     for(i = 1; i<=size(list_azero);i++)
    611622     {//rewrite the polynomials in A1 as polynomials in K[theta]
    612           tempmons = list();
    613           for (j = 1; j<=size(list_azero[i]);j++)
    614           {
    615                tempmons = tempmons + list(list_azero[i][j]);
    616           }
    617           for (j = 1 ; j<=size(tempmons);j++)
    618           {
    619                entry = leadcoef(tempmons[j]);
    620                for (k = 0; k < leadexp(tempmons[j])[2];k++)
    621                {
    622                     entry = entry*(theta-k);
    623                }
    624                tempmons[j] = entry;
    625           }
    626           list_azero[i] = sum(tempmons);
     623          tempmons = list();
     624          for (j = 1; j<=size(list_azero[i]);j++)
     625          {
     626               tempmons = tempmons + list(list_azero[i][j]);
     627          }
     628          for (j = 1 ; j<=size(tempmons);j++)
     629          {
     630               entry = leadcoef(tempmons[j]);
     631               for (k = 0; k < leadexp(tempmons[j])[2];k++)
     632               {
     633                    entry = entry*(theta-k);
     634               }
     635               tempmons[j] = entry;
     636          }
     637          list_azero[i] = sum(tempmons);
    627638     }//rewrite the polynomials in A1 as polynomials in K[theta]
    628639     //Compute all permutations of the A0-part
     
    633644     if (size(list_not_azero)!=0)
    634645     {//Compute all possibilities to permute the x's resp. the y's in the list
    635           if (list_not_azero[1] == x)
    636           {//h had a negative weighted degree
    637                shift_sign = 1;
    638                shiftvar = x;
    639           }//h had a negative weighted degree
    640           else
    641           {//h had a positive weighted degree
    642                shift_sign = -1;
    643                shiftvar = y;
    644           }//h had a positive weighted degree
    645           result = permpp(list_azero + list_not_azero);
    646           for (i = 1; i<= size(result); i++)
    647           {//adjust the a_0-parts
    648                shift = 0;
    649                for (j=1; j<=size(result[i]);j++)
    650                {
    651                     if (result[i][j]==shiftvar)
    652                     {
    653                         shift = shift + shift_sign;
    654                     }
    655                     else
    656                     {
    657                         result[i][j] = subst(result[i][j],theta,theta + shift);
    658                     }
    659                }
    660           }//adjust the a_0-parts
     646          if (list_not_azero[1] == x)
     647          {//h had a negative weighted degree
     648               shift_sign = 1;
     649               shiftvar = x;
     650          }//h had a negative weighted degree
     651          else
     652          {//h had a positive weighted degree
     653               shift_sign = -1;
     654               shiftvar = y;
     655          }//h had a positive weighted degree
     656          result = permpp(list_azero + list_not_azero);
     657          for (i = 1; i<= size(result); i++)
     658          {//adjust the a_0-parts
     659               shift = 0;
     660               for (j=1; j<=size(result[i]);j++)
     661               {
     662                    if (result[i][j]==shiftvar)
     663                    {
     664                        shift = shift + shift_sign;
     665                    }
     666                    else
     667                    {
     668                        result[i][j] = subst(result[i][j],theta,theta + shift);
     669                    }
     670               }
     671          }//adjust the a_0-parts
    661672     }//Compute all possibilities to permute the x's resp. the y's in the list
    662673     else
    663674     {//The result is just all the permutations of the a_0-part
    664           result = permpp(list_azero);
     675          result = permpp(list_azero);
    665676     }//The result is just all the permutations of the a_0 part
    666677     if (size(result)==0)
    667678     {
    668           return(result);
     679          return(result);
    669680     }
    670681     //Now we are going deeper and search for theta resp. theta + 1, substitute them by xy resp. yx and go on permuting
     
    678689     for (i = 1; i<=size(result) ; i++)
    679690     {//checking every entry of result for theta or theta +1
    680           found_theta = 0;
     691          found_theta = 0;
    681692          for(j=1;j<=size(result[i]);j++)
    682693          {
    683                if (result[i][j]==theta)
    684                {//the jth entry is theta and can be written as x*y
    685                     thetapos = j;
    686                     result[i]= insert(result[i],x,j-1);
    687                     j++;
    688                     result[i][j] = y;
    689                     found_theta = 1;
    690                     break;
    691                }//the jth entry is theta and can be written as x*y
    692                if(result[i][j] == theta +1)
    693                {
    694                     thetapos = j;
    695                     result[i] = insert(result[i],y,j-1);
    696                     j++;
    697                     result[i][j] = x;
    698                     found_theta = 1;
    699                     break;
    700                }
    701           }
    702           if (found_theta)
    703           {//One entry was theta resp. theta +1
    704                leftpart = result[i];
    705                leftpart = leftpart[1..thetapos];
    706                rightpart = result[i];
    707                rightpart = rightpart[(thetapos+1)..size(rightpart)];
    708                lparts = list(leftpart);
    709                rparts = list(rightpart);
    710                //first deal with the left part
    711                if (leftpart[thetapos] == x)
    712                {
    713                     shift_sign = 1;
    714                     shiftvar = x;
    715                }
    716                else
    717                {
    718                     shift_sign = -1;
    719                     shiftvar = y;
    720                }
    721                for (j = size(leftpart); j>1;j--)
    722                {//drip x resp. y
    723                     if (leftpart[j-1]==shiftvar)
    724                     {//commutative
    725                         j--;
    726                          continue;
    727                     }//commutative
    728                     if (deg(leftpart[j-1],intvec(-1,1,0))!=0)
    729                     {//stop here
    730                          break;
    731                     }//stop here
    732                     //Here, we can only have a a0- part
    733                     leftpart[j] = subst(leftpart[j-1],theta, theta + shift_sign);
    734                     leftpart[j-1] = shiftvar;
    735                     lparts = lparts + list(leftpart);
    736                 }//drip x resp. y
    737                 //and now deal with the right part
    738                 if (rightpart[1] == x)
    739                 {
    740                      shift_sign = 1;
    741                      shiftvar = x;
    742                 }
    743                 else
    744                 {
    745                      shift_sign = -1;
    746                      shiftvar = y;
    747                 }
    748                 for (j = 1 ; j < size(rightpart); j++)
    749                 {
    750                      if (rightpart[j+1] == shiftvar)
    751                      {
    752                           j++;
    753                           continue;
    754                      }
    755                      if (deg(rightpart[j+1],intvec(-1,1,0))!=0)
    756                      {
    757                           break;
    758                      }
    759                      rightpart[j] = subst(rightpart[j+1], theta, theta - shift_sign);
    760                      rightpart[j+1] = shiftvar;
    761                      rparts = rparts + list(rightpart);
    762                 }
    763                 //And now, we put all possibilities together
    764                 tempadd = list();
    765                 for (j = 1; j<=size(lparts); j++)
    766                 {
    767                    for (k = 1; k<=size(rparts);k++)
    768                    {
    769                               tempadd = tempadd + list(lparts[j]+rparts[k]);
    770                    }
    771                 }
    772                 tempadd = delete(tempadd,1); // The first entry is already in the list
    773                 result = result + tempadd;
    774                 continue; //We can may be not be done already with the ith entry
    775           }//One entry was theta resp. theta +1
     694               if (result[i][j]==theta)
     695               {//the jth entry is theta and can be written as x*y
     696                    thetapos = j;
     697                    result[i]= insert(result[i],x,j-1);
     698                    j++;
     699                    result[i][j] = y;
     700                    found_theta = 1;
     701                    break;
     702               }//the jth entry is theta and can be written as x*y
     703               if(result[i][j] == theta +1)
     704               {
     705                    thetapos = j;
     706                    result[i] = insert(result[i],y,j-1);
     707                    j++;
     708                    result[i][j] = x;
     709                    found_theta = 1;
     710                    break;
     711               }
     712          }
     713          if (found_theta)
     714          {//One entry was theta resp. theta +1
     715               leftpart = result[i];
     716               leftpart = leftpart[1..thetapos];
     717               rightpart = result[i];
     718               rightpart = rightpart[(thetapos+1)..size(rightpart)];
     719               lparts = list(leftpart);
     720               rparts = list(rightpart);
     721               //first deal with the left part
     722               if (leftpart[thetapos] == x)
     723               {
     724                    shift_sign = 1;
     725                    shiftvar = x;
     726               }
     727               else
     728               {
     729                    shift_sign = -1;
     730                    shiftvar = y;
     731               }
     732               for (j = size(leftpart); j>1;j--)
     733               {//drip x resp. y
     734                    if (leftpart[j-1]==shiftvar)
     735                    {//commutative
     736                        j--;
     737                         continue;
     738                    }//commutative
     739                    if (deg(leftpart[j-1],intvec(-1,1,0))!=0)
     740                    {//stop here
     741                         break;
     742                    }//stop here
     743                    //Here, we can only have a a0- part
     744                    leftpart[j] = subst(leftpart[j-1],theta, theta + shift_sign);
     745                    leftpart[j-1] = shiftvar;
     746                    lparts = lparts + list(leftpart);
     747                }//drip x resp. y
     748                //and now deal with the right part
     749                if (rightpart[1] == x)
     750                {
     751                     shift_sign = 1;
     752                     shiftvar = x;
     753                }
     754                else
     755                {
     756                     shift_sign = -1;
     757                     shiftvar = y;
     758                }
     759                for (j = 1 ; j < size(rightpart); j++)
     760                {
     761                     if (rightpart[j+1] == shiftvar)
     762                     {
     763                          j++;
     764                          continue;
     765                     }
     766                     if (deg(rightpart[j+1],intvec(-1,1,0))!=0)
     767                     {
     768                          break;
     769                     }
     770                     rightpart[j] = subst(rightpart[j+1], theta, theta - shift_sign);
     771                     rightpart[j+1] = shiftvar;
     772                     rparts = rparts + list(rightpart);
     773                }
     774                //And now, we put all possibilities together
     775                tempadd = list();
     776                for (j = 1; j<=size(lparts); j++)
     777                {
     778                        for (k = 1; k<=size(rparts);k++)
     779                   {
     780                              tempadd = tempadd + list(lparts[j]+rparts[k]);
     781                   }
     782                }
     783                tempadd = delete(tempadd,1); // The first entry is already in the list
     784                result = result + tempadd;
     785                continue; //We can may be not be done already with the ith entry
     786          }//One entry was theta resp. theta +1
    776787     }//checking every entry of result for theta or theta +1
    777788     //map back to the basering
     
    781792     for (i=1; i<=size(result);i++)
    782793     {//adding the K factor
    783           result[i] = k_factor + result[i];
     794          result[i] = k_factor + result[i];
    784795     }//adding the k-factor
    785796     result = delete_dublicates_noteval(result);
     
    805816     if (size(l)==0)
    806817     {
    807           return(list());
     818          return(list());
    808819     }
    809820     if (size(l)==1)
    810821     {
    811           return(list(l));
     822          return(list(l));
    812823     }
    813824     for (i = 1; i<=size(l); i++ )
    814825     {
    815           tempresult = perm(delete(l,i));
    816           for (j = 1; j<=size(tempresult);j++)
    817           {
    818                tempresult[j] = list(l[i])+tempresult[j];
    819           }
    820           result = result+tempresult;
     826          tempresult = perm(delete(l,i));
     827          for (j = 1; j<=size(tempresult);j++)
     828          {
     829               tempresult[j] = list(l[i])+tempresult[j];
     830          }
     831          result = result+tempresult;
    821832     }
    822833     return(result);
     
    836847     if (size(l)==0)
    837848     {
    838           return(list());
     849          return(list());
    839850     }
    840851     if (size(l)==1)
    841852     {
    842           return(list(l));
     853          return(list(l));
    843854     }
    844855     for (i = 1; i<=size(l);i++)
    845856     {//Filling the list with unique entries
    846           double_entry = 0;
    847           for (j = 1; j<=size(l_without_double);j++)
    848           {
    849                if (l_without_double[j] == l[i])
    850                {
    851                     double_entry = 1;
    852                     break;
    853                }
    854           }
    855           if (!double_entry)
    856           {
    857                l_without_double = l_without_double + list(l[i]);
    858                l_without_double_pos = l_without_double_pos + list(i);
    859           }
     857          double_entry = 0;
     858          for (j = 1; j<=size(l_without_double);j++)
     859          {
     860               if (l_without_double[j] == l[i])
     861               {
     862                    double_entry = 1;
     863                    break;
     864               }
     865          }
     866          if (!double_entry)
     867          {
     868               l_without_double = l_without_double + list(l[i]);
     869               l_without_double_pos = l_without_double_pos + list(i);
     870          }
    860871     }//Filling the list with unique entries
    861872     for (i = 1; i<=size(l_without_double); i++ )
    862873     {
    863           tempresult = permpp(delete(l,l_without_double_pos[i]));
    864           for (j = 1; j<=size(tempresult);j++)
    865           {
    866                tempresult[j] = list(l_without_double[i])+tempresult[j];
    867           }
    868           result = result+tempresult;
     874          tempresult = permpp(delete(l,l_without_double_pos[i]));
     875          for (j = 1; j<=size(tempresult);j++)
     876          {
     877               tempresult[j] = list(l_without_double[i])+tempresult[j];
     878          }
     879          result = result+tempresult;
    869880     }
    870881     return(result);
     
    878889//the procedure sfacwa, where the ring must contain the
    879890//variables in a certain order.
    880 proc facfirstwa(poly h)
    881 {//proc facfirstwa
     891proc facFirstWeyl(poly h)
     892"USAGE: facFirstWeyl(h); h a poly
     893RETURN: list
     894PURPOSE: compute all factorizations of a polynomial in the first Weyl algebra
     895THEORY: Implements the new algorithm by A. Heinle and V. Levandovskyy, see the thesis of A. Heinle
     896ASSUME: basering in the first Weyl algebra
     897EXAMPLE: example facFirstWeyl; shows examples
     898SEE ALSO: facSubWeyl, testNCfac
     899"{//proc facFirstWeyl
    882900      //Redefine the ring in my standard form
    883901     if (!isWeyl())
    884902     {//Our basering is not the Weyl algebra
    885           return(list());
     903          return(list());
    886904     }//Our basering is not the Weyl algebra
    887905     if(nvars(basering)!=2)
    888906     {//Our basering is the Weyl algebra, but not the first
    889           return(list());
     907          return(list());
    890908     }//Our basering is the Weyl algebra, but not the first
    891909     if (ringlist(basering)[6][1,2] == -1) //manual of ringlist will tell you why
    892910     {
    893           def r = basering;
    894           ring tempRing = ringlist(r)[1][1],(x,y),Ws(-1,1);
    895           def tempRingnc = nc_algebra(1,1);
    896           setring(tempRingnc);
    897           map transf = r, var(2), var(1);
    898           list result = sfacwa(transf(h));
    899           setring(r);
    900           map transfback = tempRingnc, var(2),var(1);
    901           return(transfback(result));
     911          def r = basering;
     912          ring tempRing = ringlist(r)[1][1],(x,y),Ws(-1,1);
     913          def tempRingnc = nc_algebra(1,1);
     914          setring(tempRingnc);
     915          map transf = r, var(2), var(1);
     916          list result = sfacwa(transf(h));
     917          setring(r);
     918          map transfback = tempRingnc, var(2),var(1);
     919          return(transfback(result));
    902920     }
    903921     else { return(sfacwa(h));}
    904 }//proc facfirstwa
     922}//proc facFirstWeyl
    905923example
    906924{
     
    909927     def r = nc_algebra(1,1);
    910928     setring(r);
    911      poly h = (x^2*y^2+x)*(y-1);
    912      facfirstwa(h);
     929     poly h = (x^2*y^2+x)*(x+1);
     930     facFirstWeyl(h);
    913931}
    914932
    915933//This is the main program
    916934static proc sfacwa(poly h)
    917 "USAGE: facwa(h); h is a polynomial in the first Weyl algebra
     935"USAGE: sfacwa(h); h is a polynomial in the first Weyl algebra
    918936RETURN: list
    919937PURPOSE: Computes a factorization of a polynomial h in the first Weyl algebra
    920 THEORY: @code{homogfac} returns a list with some factorizations of the given polynomial. The possibilities of the factorization of the highest homogeneous part and those of the lowest will be merged.
    921 EXAMPLE: example facwa; shows examples
     938THEORY: @code{sfacwa} returns a list with some factorizations of the given
     939        polynomial. The possibilities of the factorization of the highest
     940        homogeneous part and those of the lowest will be merged. If during this
     941        procedure a factorization of the polynomial occurs, it will be added to
     942        the output list. For a more detailled description visit
     943        @url{http://math.rwth-aachen.de/\~Albert.Heinle}
    922944SEE ALSO: homogfac_all, homogfac
    923945"{//proc sfacwa
    924946     if(homogwithorder(h,intvec(-1,1)))
    925947     {
    926           return(homogfac_all(h));
     948          return(homogfac_all(h));
    927949     }
    928950     def r = basering;
     
    944966     for (i = 1; i<=size(f1);i++)
    945967     {//Merge all combinations
    946           for (j = 1; j<=size(f2); j++)
    947           {
    948                M = M + mergence(f1[i],f2[j],limits);
    949           }
     968          for (j = 1; j<=size(f2); j++)
     969          {
     970               M = M + mergence(f1[i],f2[j],limits);
     971          }
    950972     }//Merge all combinations
    951973     for (i = 1 ; i<= size(M); i++)
    952974     {//filter valid combinations
    953           if (product(M[i]) == h)
    954           {//We have one factorization
    955                result = result + list(M[i]);
    956                M = delete(M,i);
    957                continue;
    958           }//We have one factorization
    959           else
    960           {
    961                if (deg(h,intvec(-1,1))<=deg(h-product(M[i]),intvec(-1,1)))
    962                {
    963                     M = delete(M,i);
    964                     continue;
    965                }
    966                if (deg(h,intvec(1,-1))<=deg(h-product(M[i]),intvec(1,-1)))
    967                {
    968                     M = delete(M,i);
    969                     continue;
    970                }
    971           }
     975          if (product(M[i]) == h)
     976          {//We have one factorization
     977               result = result + list(M[i]);
     978               M = delete(M,i);
     979               continue;
     980          }//We have one factorization
     981          else
     982          {
     983               if (deg(h,intvec(-1,1))<=deg(h-product(M[i]),intvec(-1,1)))
     984               {
     985                    M = delete(M,i);
     986                    continue;
     987               }
     988               if (deg(h,intvec(1,-1))<=deg(h-product(M[i]),intvec(1,-1)))
     989               {
     990                    M = delete(M,i);
     991                    continue;
     992               }
     993          }
    972994     }//filter valid combinations
    973995     M = delete_dublicates_eval(M);
    974996     while(size(M)>0)
    975997     {//work on the elements of M
    976           hatM = list();
    977           for(i = 1; i<=size(M); i++)
    978           {//iterate over all elements of M
    979                hath = h-product(M[i]);
    980                temp = list();
    981                //First check for common inhomogeneous factors between hath and h
    982                if (involution(NF(involution(hath,invo), std(involution(ideal(M[i][1]),invo))),invo)==0)
    983                {//hath and h have a common factor on the left
    984                     j = 1;
    985                     f1 = M[i];
    986                     if (j+1<=size(f1))
    987                     {//Checking for more than one common factor
    988                         while(involution(NF(involution(hath,invo),std(involution(ideal(product(f1[1..(j+1)])),invo))),invo)==0)
    989                         {
    990                               if (j+1<size(f1))
    991                               {
    992                                    j++;
    993                               }
    994                               else
    995                               {
    996                                    break;
    997                               }
    998                         }
    999                     }//Checking for more than one common factor
    1000                     f2 = list(f1[1..j])+list(involution(lift(involution(product(f1[1..j]),invo),involution(hath,invo))[1,1],invo));
    1001                     temp = temp + merge_cf(f2,f1,limits);
    1002                }//hath and h have a common factor on the left
    1003                if (reduce(hath, std(ideal(M[i][size(M[i])])))==0)
    1004                {//hath and h have a common factor on the right
    1005                     j = size(M[i]);
    1006                     f1 = M[i];
    1007                     if (j-1>0)
    1008                     {//Checking for more than one factor
    1009                         while(reduce(hath,std(ideal(product(f1[(j-1)..size(f1)]))))==0)
    1010                         {
    1011                               if (j-1>1)
    1012                               {
    1013                                    j--;
    1014                               }
    1015                               else
    1016                               {
    1017                                    break;
    1018                               }
    1019                         }
    1020                     }//Checking for more than one factor
    1021                     f2 = list(lift(product(f1[j..size(f1)]),hath)[1,1])+list(f1[j..size(f1)]);
    1022                     temp = temp + merge_cf(f2,M[i],limits);
    1023                }//hath and h have a common factor on the right
    1024                //and now the homogeneous
    1025                maxh = jet(hath,deg(hath,intvec(-1,1)),intvec(-1,1))-jet(hath,deg(hath,intvec(-1,1))-1,intvec(-1,1));
    1026                minh = jet(hath,deg(hath,intvec(1,-1)),intvec(1,-1))-jet(hath,deg(hath,intvec(1,-1))-1,intvec(1,-1));
    1027                f1 = homogfac_all(maxh);
    1028                f2 = homogfac_all(minh);
    1029                for (j = 1; j<=size(f1);j++)
    1030                {
    1031                     for (k=1; k<=size(f2);k++)
    1032                     {
    1033                         homogtemp = mergence(f1[j],f2[k],limits);
    1034                     }
    1035                }
    1036                for (j = 1; j<= size(homogtemp); j++)
    1037                {
    1038                     temp = temp + mergence(homogtemp[j],M[i],limits);
    1039                }
    1040                for (j = 1; j<=size(temp); j++)
    1041                {//filtering invalid entries in temp
    1042                     if(product(temp[j])==h)
    1043                     {//This is already a result
    1044                         result = result + list(temp[j]);
    1045                         temp = delete(temp,j);
    1046                         continue;
    1047                     }//This is already a result
    1048                     if (deg(hath,intvec(-1,1))<=deg(hath-product(temp[j]),intvec(-1,1)))
    1049                     {
    1050                         temp = delete(temp,j);
    1051                         continue;
    1052                     }
    1053                }//filtering invalid entries in temp
    1054                hatM = hatM + temp;
    1055           }//iterate over all elements of M
    1056           M = hatM;
    1057           for (i = 1; i<=size(M);i++)
    1058           {//checking for complete factorizations
    1059                if (h == product(M[i]))
    1060                {
    1061                     result = result + list(M[i]);
    1062                     M = delete(M,i);
    1063                     continue;
    1064                }
    1065           }//checking for complete factorizations
    1066           M = delete_dublicates_eval(M);
     998          hatM = list();
     999          for(i = 1; i<=size(M); i++)
     1000          {//iterate over all elements of M
     1001               hath = h-product(M[i]);
     1002               temp = list();
     1003               //First check for common inhomogeneous factors between hath and h
     1004               if (involution(NF(involution(hath,invo), std(involution(ideal(M[i][1]),invo))),invo)==0)
     1005               {//hath and h have a common factor on the left
     1006                    j = 1;
     1007                    f1 = M[i];
     1008                    if (j+1<=size(f1))
     1009                    {//Checking for more than one common factor
     1010                        while(involution(NF(involution(hath,invo),std(involution(ideal(product(f1[1..(j+1)])),invo))),invo)==0)
     1011                        {
     1012                              if (j+1<size(f1))
     1013                              {
     1014                                   j++;
     1015                              }
     1016                              else
     1017                              {
     1018                                   break;
     1019                              }
     1020                        }
     1021                    }//Checking for more than one common factor
     1022                    f2 = list(f1[1..j])+list(involution(lift(involution(product(f1[1..j]),invo),involution(hath,invo))[1,1],invo));
     1023                    temp = temp + merge_cf(f2,f1,limits);
     1024               }//hath and h have a common factor on the left
     1025               if (reduce(hath, std(ideal(M[i][size(M[i])])))==0)
     1026               {//hath and h have a common factor on the right
     1027                    j = size(M[i]);
     1028                    f1 = M[i];
     1029                    if (j-1>0)
     1030                    {//Checking for more than one factor
     1031                        while(reduce(hath,std(ideal(product(f1[(j-1)..size(f1)]))))==0)
     1032                        {
     1033                              if (j-1>1)
     1034                              {
     1035                                   j--;
     1036                              }
     1037                              else
     1038                              {
     1039                                   break;
     1040                              }
     1041                        }
     1042                    }//Checking for more than one factor
     1043                         f2 = list(lift(product(f1[j..size(f1)]),hath)[1,1])+list(f1[j..size(f1)]);
     1044                    temp = temp + merge_cf(f2,M[i],limits);
     1045               }//hath and h have a common factor on the right
     1046               //and now the homogeneous
     1047               maxh = jet(hath,deg(hath,intvec(-1,1)),intvec(-1,1))-jet(hath,deg(hath,intvec(-1,1))-1,intvec(-1,1));
     1048               minh = jet(hath,deg(hath,intvec(1,-1)),intvec(1,-1))-jet(hath,deg(hath,intvec(1,-1))-1,intvec(1,-1));
     1049               f1 = homogfac_all(maxh);
     1050               f2 = homogfac_all(minh);
     1051               for (j = 1; j<=size(f1);j++)
     1052               {
     1053                    for (k=1; k<=size(f2);k++)
     1054                    {
     1055                        homogtemp = mergence(f1[j],f2[k],limits);
     1056                    }
     1057               }
     1058               for (j = 1; j<= size(homogtemp); j++)
     1059               {
     1060                    temp = temp + mergence(homogtemp[j],M[i],limits);
     1061               }
     1062               for (j = 1; j<=size(temp); j++)
     1063               {//filtering invalid entries in temp
     1064                    if(product(temp[j])==h)
     1065                    {//This is already a result
     1066                        result = result + list(temp[j]);
     1067                        temp = delete(temp,j);
     1068                        continue;
     1069                    }//This is already a result
     1070                    if (deg(hath,intvec(-1,1))<=deg(hath-product(temp[j]),intvec(-1,1)))
     1071                    {
     1072                        temp = delete(temp,j);
     1073                        continue;
     1074                    }
     1075               }//filtering invalid entries in temp
     1076               hatM = hatM + temp;
     1077          }//iterate over all elements of M
     1078          M = hatM;
     1079          for (i = 1; i<=size(M);i++)
     1080          {//checking for complete factorizations
     1081               if (h == product(M[i]))
     1082               {
     1083                    result = result + list(M[i]);
     1084                    M = delete(M,i);
     1085                    continue;
     1086               }
     1087          }//checking for complete factorizations
     1088          M = delete_dublicates_eval(M);
    10671089     }//work on the elements of M
    10681090     //In the case, that there is none, write a constant factor before the factor of interest.
    10691091     for (i = 1 ; i<=size(result);i++)
    10701092     {//add a constant factor
    1071           if (deg(result[i][1],intvec(1,1))!=0)
    1072           {
    1073                result[i] = insert(result[i],1);
    1074           }
     1093          if (deg(result[i][1],intvec(1,1))!=0)
     1094          {
     1095               result[i] = insert(result[i],1);
     1096          }
    10751097     }//add a constant factor
    10761098     result = delete_dublicates_noteval(result);
     
    10871109     if(size(weights) != nvars(basering))
    10881110     {//The user does not know how many variables the current ring has
    1089           return(0);
     1111          return(0);
    10901112     }//The user does not know how many variables the current ring has
    10911113     int i;
     
    10931115     for (i = 1; i<=size(h);i++)
    10941116     {
    1095           if (deg(h[i],weights)!=dofp)
    1096           {
    1097                return(0);
    1098           }
     1117          if (deg(h[i],weights)!=dofp)
     1118          {
     1119               return(0);
     1120          }
    10991121     }
    11001122     return(1);
     
    11131135// and if not, a list with the products subtracted by p
    11141136// will be returned
    1115 proc testfac(list l, list #)
    1116 "USAGE:  - testfac(l); l is a list, that contains lists with entries in the first Weyl algebra or
    1117          @*- testfac(l,p); l is a list, that contains lists with entries in the first Weyl algebra and p is a polynomial in the first Weyl algebra
    1118 RETURN: list@*
    1119 PURPOSE: Checks whether a list of factorizations contains factorizations of the same element in the first Weyl algebra@*
    1120 THEORY: @code{testfac} multiplies out each factorization and checks whether each factorization was a factorization of the same element.
    1121 @* - if there is only a list given, the output will be the empty list, if it does not contain factorizations of the same element. Otherwise it will multiply out all the factorizations and return a list of the length if the given one, in which each entry is the factorized polynomial
    1122 @* - if there is a polynomial in the second argument, then the procedure checks whether the given list contains factorizations of this polynomial. If it does, then the output will be a list with the length of the given one and each entry contains this polynomial. If it does not, then the procedure returns a list, in which each entry is the given polynomial subtracted by the multiplied candidate for its factorization.
    1123 EXAMPLE: example testfac; shows examples
    1124 SEE ALSO: homogfac_all, homogfac, facwa
     1137proc testNCfac(list l, list #)
     1138"USAGE: testNCfac(l[,p]); l is a list, p is an optional poly
     1139RETURN: list
     1140ASSUME: basering is the first Weyl algebra, the entries of l are polynomials
     1141PURPOSE: Checks whether a list of factorizations contains factorizations of
     1142         the same element in the first Weyl algebra
     1143THEORY: @code{testNCfac} multiplies out each factorization and checks whether
     1144        each factorization was a factorization of the same element.
     1145@* - if there is only a list given, the output will be the empty list, if it
     1146     does not contain factorizations of the same element. Otherwise it will
     1147     multiply out all the factorizations and return a list of the length if the
     1148     given one, in which each entry is the factorized polynomial
     1149@* - if there is a polynomial in the second argument, then the procedure checks
     1150     whether the given list contains factorizations of this polynomial. If it
     1151     does, then the output will be a list with the length of the given one and
     1152     each entry contains this polynomial. If it does not, then the procedure
     1153     returns a list, in which each entry is the given polynomial subtracted
     1154     by the multiplied candidate for its factorization.
     1155EXAMPLE: example testNCfac; shows examples
     1156SEE ALSO: facFirstWeyl, facSubWeyl
    11251157"{//proc testfac
    11261158     if (size(l)==0)
    11271159     {//The empty list is given
    1128           return(list());
     1160          return(list());
    11291161     }//The empty list is given
    11301162     if (size(#)>1)
    11311163     {//We want max. one optional argument
    1132           return(list());
     1164          return(list());
    11331165     }//We want max. one optional argument
    11341166     list result;
     
    11361168     if (size(#)==0)
    11371169     {//No optional argument is given
    1138           int valid = 1;
    1139           for (i = size(l);i>=1;i--)
    1140           {//iterate over the elements of the given list
    1141                if (size(result)>0)
    1142                {
    1143                     if (product(l[i])!=result[size(l)-i])
    1144                     {
    1145                         valid = 0;
    1146                         break;
    1147                     }
    1148                }
    1149                result = insert(result, product(l[i]));
    1150           }//iterate over the elements of the given list
    1151           if (valid)
    1152           {
    1153                return(result);
    1154           }
    1155           return(list());
     1170          int valid = 1;
     1171          for (i = size(l);i>=1;i--)
     1172          {//iterate over the elements of the given list
     1173               if (size(result)>0)
     1174               {
     1175                    if (product(l[i])!=result[size(l)-i])
     1176                    {
     1177                        valid = 0;
     1178                        break;
     1179                    }
     1180               }
     1181               result = insert(result, product(l[i]));
     1182          }//iterate over the elements of the given list
     1183          if (valid)
     1184          {
     1185               return(result);
     1186          }
     1187          return(list());
    11561188     }//No optional argument is given
    11571189     else
    11581190     {
    1159           int valid = 1;
    1160           for (i = size(l);i>=1;i--)
    1161           {//iterate over the elements of the given list
    1162                if (product(l[i])!=#[1])
    1163                {
    1164                     valid = 0;
    1165                }
    1166                result = insert(result, product(l[i]));
    1167           }//iterate over the elements of the given list
    1168           if (valid)
    1169           {
    1170                return(result);
    1171           }
    1172           for (i = 1 ; i<= size(result);i++)
    1173           {//subtract p from every entry in result
    1174                result[i] = result[i] - #[1];
    1175           }//subtract p from every entry in result
    1176           return(result);
     1191          int valid = 1;
     1192          for (i = size(l);i>=1;i--)
     1193          {//iterate over the elements of the given list
     1194               if (product(l[i])!=#[1])
     1195               {
     1196                    valid = 0;
     1197               }
     1198               result = insert(result, product(l[i]));
     1199          }//iterate over the elements of the given list
     1200          if (valid)
     1201          {
     1202               return(result);
     1203          }
     1204          for (i = 1 ; i<= size(result);i++)
     1205          {//subtract p from every entry in result
     1206               result[i] = result[i] - #[1];
     1207          }//subtract p from every entry in result
     1208          return(result);
    11771209     }
    11781210}//proc testfac
     
    11801212{
    11811213     "EXAMPLE:";echo=2;
    1182      ring R = 0,(x,y),Ws(-1,1);
    1183      def r = nc_algebra(1,1);
    1184      setring(r);
    1185      poly h = (x^2*y^2+1)*(x^4);
    1186      def t1 = facfirstwa(h);
     1214     ring r = 0,(x,y),dp;
     1215     def R = nc_algebra(1,1);
     1216     setring R;
     1217     poly h = (x^2*y^2+1)*(x^2);
     1218     def t1 = facFirstWeyl(h);
    11871219     //fist a correct list
    1188      testfac(t1);
     1220     testNCfac(t1);
    11891221     //now a correct list with the factorized polynomial
    1190      testfac(t1,h);
     1222     testNCfac(t1,h);
    11911223     //now we put in an incorrect list without a polynomial
    1192      t1[5][3] = y;
    1193      testfac(t1);
     1224     t1[3][3] = y;
     1225     testNCfac(t1);
    11941226     //and last but not least we take h as additional input
    1195      testfac(t1,h);
     1227     testNCfac(t1,h);
    11961228}
    11971229//==================================================
    1198 //Procedure facsubwa:
     1230//Procedure facSubWeyl:
    11991231//This procedure serves the purpose to compute a
    12001232//factorization of a given polynomial in a ring, whose subring
    12011233//is the first Weyl algebra. The polynomial must only contain
    12021234//the two arguments, which are also given by the user.
    1203 proc facsubwa(poly h, X, D)
    1204 "USAGE:  facsubwa(h,x,y) - h is a polynomial in the first Weyl algebra, X and D are the variables in h, which have the commutation rule DX = XD + 1@*
    1205 RETURN: list@*
    1206 PURPOSE: Given a polynomial in a Ring, which contains the first Weyl algebra as a subring. Furthermore let this polynomial be in the first Weyl algebra with the variables X and D with DX = XD +1.@*
    1207 THEORY: @code{facsubwa} computes some factorizations of the given polynomial, which is a polynomial in the variables X and D@*
    1208 EXAMPLE: example facsubwa; shows examples
    1209 SEE ALSO: facfirstwa
    1210 "{//proc facsubwa
    1211      if(!isWeyl())
    1212      {
    1213           return(list());
    1214      }
     1235
     1236proc facSubWeyl(poly h, X, D)
     1237"USAGE:  facSubWeyl(h,x,y); h, X, D polynomials
     1238RETURN: list
     1239ASSUME: X and D are variables of a basering, which satisfy DX = XD +1.
     1240@* That is,  they generate the copy of the first Weyl algebra in a basering.
     1241@* Moreover, h is a polynomial in X and D only.
     1242PURPOSE: compute factorizations of the polynomial, which depends on X and D.
     1243EXAMPLE: example facSubWeyl; shows examples
     1244SEE ALSO: facFirstWeyl, testNCfac
     1245"{
     1246//proc facSubWeyl
     1247// basering can be anything having a Weyl algebra as subalgebra
    12151248     def r = basering;
    1216      //We begin to check the input
    1217      if (size(X)>1)
    1218      {
    1219           return(list());
    1220      }
    1221      if (size(D)>1)
    1222      {
    1223           return(list());
    1224      }
     1249     //We begin to check the input for assumptions
     1250     // which are: X,D are vars of the basering,
     1251     if ( (isVar(X)!=1) || (isVar(D)!=1) || (size(X)>1) || (size(D)>1) ||
     1252          (leadcoef(X) != number(1)) || (leadcoef(D) != number(1)) )
     1253     {
     1254       ERROR("expected pure variables as generators of a subalgebra");
     1255     }
     1256     // Weyl algebra:
     1257     poly w = D*X-X*D-1; // [D,X]=1
     1258     poly u = D*X-X*D+1; // [X,D]=1
     1259     if (u*w!=0)
     1260     {
     1261       // that is no combination gives Weyl
     1262       ERROR("2nd and 3rd argument do not generate a Weyl algebra");
     1263     }
     1264     // one of two is correct
     1265     if (u==0)
     1266     {
     1267       // hence w != 0, swap X<->D
     1268       w = D; D=X; X=w;
     1269     }
     1270     // DONE with assumptions     
     1271     //Input successfully checked
    12251272     intvec lexpofX = leadexp(X);
    12261273     intvec lexpofD = leadexp(D);
    1227      if (sum(lexpofX)!=1)
    1228      {
    1229           return(list());
    1230      }
    1231      if (sum(lexpofD)!=1)
    1232      {
    1233           return(list());
    1234      }
    12351274     int varnumX=1;
    12361275     int varnumD=1;
    12371276     while(lexpofX[varnumX] != 1)
    12381277     {
    1239           varnumX++;
     1278          varnumX++;
    12401279     }
    12411280     while(lexpofD[varnumD] != 1)
    12421281     {
    1243           varnumD++;
    1244      }
    1245      lexpofX = lexpofX +1;
    1246      lexpofX[varnumX] = 0;
    1247      lexpofX[varnumD] = 0;
    1248      if (deg(h,lexpofX)!=0)
    1249      {
    1250           return(list());
    1251      }
    1252      //Input successfully checked
     1282          varnumD++;
     1283     }
     1284/* VL: it's not clear what you do here!!! I comment out */
     1285     //     lexpofX = lexpofX +1;
     1286     //      lexpofX[varnumX] = 0;
     1287     //      lexpofX[varnumD] = 0;
     1288     //      if (deg(h,lexpofX)!=0)
     1289     //      {
     1290     //           return(list());
     1291     //      }
     1292/* VL : to add printlevel stuff */
     1293
    12531294     ring firstweyl = 0,(var(varnumX),var(varnumD)),dp;
    12541295     def Firstweyl = nc_algebra(1,1);
    12551296     setring Firstweyl;
    12561297     poly h = imap(r,h);
    1257      list result = facfirstwa(h);
     1298     list result = facFirstWeyl(h);
    12581299     setring r;
    12591300     list result = imap(Firstweyl,result);
    12601301     return(result);
    1261 }//proc facsubwa
     1302}//proc facSubWeyl
    12621303example
    12631304{
    12641305     "EXAMPLE:";echo=2;
    1265      ring R = 0,(x,y),Ws(-1,1);
    1266      def r = nc_algebra(1,1);
    1267      setring(r);
    1268      poly h = (x^2*y^2+x)*x;
    1269      facsubwa(h,x,y);
     1306     ring r = 0,(x,y,z),dp;
     1307     matrix D[3][3]; D[1,3]=-1;
     1308     def R = nc_algebra(1,D); // x,z generate Weyl subalgebra
     1309     setring R;
     1310     poly h = (x^2*z^2+x)*x;
     1311     list fact1 = facSubWeyl(h,x,z); fact1[2];
     1312     poly test1 = fact1[2][2]; // for testing
     1313     // compare with facFirstWeyl:
     1314     ring s = 0,(z,x),dp;
     1315     def S = nc_algebra(1,1); setring S;
     1316     poly h = (x^2*z^2+x)*x;
     1317     list fact2 = facFirstWeyl(h); fact2[2];
     1318     poly test2 = fact2[2][2]; // for testing
     1319     map F = R,x,z;
     1320     F(test1) - test2; // they are identical
    12701321}
    12711322//==================================================
Note: See TracChangeset for help on using the changeset viewer.