Changeset cfb500 in git


Ignore:
Timestamp:
May 20, 2011, 1:07:49 PM (13 years ago)
Author:
Frank Seelisch <seelisch@…>
Branches:
(u'fieker-DuVal', '117eb8c30fc9e991c4decca4832b1d19036c4c65')(u'spielwiese', '45e267b9942dec9429fe249ce3e5e44ab80a6a3a')
Children:
502f5f6f3100e91a41cf8b3762727582bcdcc2fa
Parents:
09ed7738e7d4dd96c01950a6f0736d8581f333d2
git-author:
Frank Seelisch <seelisch@mathematik.uni-kl.de>2011-05-20 13:07:49+02:00
git-committer:
Mohamed Barakat <mohamed.barakat@rwth-aachen.de>2011-11-09 12:36:13+01:00
Message:
ending tests for alg field ext impl (all tests ok)
Location:
libpolys
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • libpolys/polys/ext_fields/algext.cc

    r09ed77 rcfb500  
    437437}
    438438
     439/* IMPORTANT NOTE: Since an algebraic field extension is again a field,
     440                   the gcd of two elements is not very interesting. (It
     441                   is actually any unit in the field, i.e., any non-
     442                   zero element.) Note that the below method does not operate
     443                   in this strong sense but rather computes the gcd of
     444                   two given elements in the underlying polynomial ring. */
    439445number naGcd(number a, number b, const coeffs cf)
    440446{
     
    455461  p_Delete(&theGcd, naRing);
    456462  p_Delete(&mFactor, naRing);
    457   /* printf("naInvers\n");
    458      p_Write((poly)a, naRing);
    459      p_Write(aFactor, naRing); */
    460463  return (number)(aFactor);
    461464}
  • libpolys/polys/monomials/p_polys.cc

    r09ed77 rcfb500  
    14451445poly p_PolyDiv(poly &p, poly divisor, BOOLEAN needResult, ring r)
    14461446{
    1447 //printf("p_PolyDiv:\n");
    14481447  assume(divisor != NULL);
    14491448  if (p == NULL) return NULL;
     
    14521451  number divisorLC = p_GetCoeff(divisor, r);
    14531452  int divisorLE = p_GetExp(divisor, 1, r);
    1454 //p_Write(p, r); p_Write(divisor, r);
    14551453  while ((p != NULL) && (p_Deg(p, r) >= p_Deg(divisor, r)))
    14561454  {
     
    14621460    p_SetExp(t, 1, e, r);
    14631461    p_Setm(t, r);
    1464 //printf("t\n");
    1465 //p_Write(t, r);
    14661462    if (needResult) result = p_Add_q(result, p_Copy(t, r), r);
    14671463    p = p_Add_q(p, p_Neg(p_Mult_q(t, p_Copy(divisor, r), r), r), r);
    1468 //printf("p\n");
    1469 //p_Write(p, r);
    14701464  }
    14711465  return result;
     
    14831477  number lc = p_GetCoeff(p, r);
    14841478  if (n_IsOne(lc, r->cf)) return;
    1485 
    1486   if( pNext(p) == NULL )
    1487   {
    1488      p_SetCoeff(p, n_Init(1, r->cf), r); // will delete the leading coeff 1st!
    1489      return;
    1490   }
    1491    
    1492   // TODO: consider unsing 1xn_Inverse + n_Mult instead of multiple n_Div
    1493   p_SetCoeff0(p, n_Init(1, r->cf), r); // no coeff destruction!
    1494    
     1479  number lcInverse = n_Invers(lc, r->cf);
     1480  number n = n_Init(1, r->cf);
     1481  p_SetCoeff(p, n, r);   // destroys old leading coefficient!
    14951482  p = pIter(p);
    1496    
    14971483  while (p != NULL)
    14981484  {
    1499     number c = p_GetCoeff(p, r);
    1500     number n = n_Div(c, lc, r->cf);
    1501     n_Delete(&c, r->cf);
    1502     p_SetCoeff(p, n, r);
     1485    number n = n_Mult(p_GetCoeff(p, r), lcInverse, r->cf);
     1486    p_SetCoeff(p, n, r);   // destroys old leading coefficient!
    15031487    p = pIter(p);
    15041488  }
    1505   n_Delete(&lc, r->cf);
     1489  n_Delete(&lcInverse, r->cf);
    15061490  p = pp;
    15071491}
     
    15521536    qFactor = NULL;
    15531537    pFactor = p_ISet(1, r);
    1554     number n = p_GetCoeff(pFactor, r);
    1555 //printf("p_ExtGcdHelper0:\n");
    1556 //p_Write(p, r);
    15571538    p_SetCoeff(pFactor, n_Invers(p_GetCoeff(p, r), r->cf), r);
    1558     n_Delete(&n, r->cf);
    15591539    p_Monic(p, r);
    15601540    return p;
     
    15621542  else
    15631543  {
    1564 //printf("p_ExtGcdHelper1:\n");
    1565 //p_Write(p, r); p_Write(q, r);
    15661544    poly pDivQ = p_PolyDiv(p, q, TRUE, r);
    15671545    poly ppFactor = NULL; poly qqFactor = NULL;
    1568 //printf("p_ExtGcdHelper2:\n");
    1569 //p_Write(p, r); p_Write(ppFactor, r);
    1570 //p_Write(q, r); p_Write(qqFactor, r);
    15711546    poly theGcd = p_ExtGcdHelper(q, qqFactor, p, ppFactor, r);
    1572 //printf("p_ExtGcdHelper3:\n");
    1573 //p_Write(q, r); p_Write(qqFactor, r);
    1574 //p_Write(p, r); p_Write(ppFactor, r);
    1575 //p_Write(theGcd, r);
    15761547    pFactor = ppFactor;
    15771548    qFactor = p_Add_q(qqFactor,
    15781549                      p_Neg(p_Mult_q(pDivQ, p_Copy(ppFactor, r), r), r),
    15791550                      r);
    1580 //printf("p_ExtGcdHelper4:\n");
    1581 //p_Write(pFactor, r); p_Write(qFactor, r);
    15821551    return theGcd;
    15831552  }
     
    15941563poly p_ExtGcd(poly p, poly &pFactor, poly q, poly &qFactor, ring r)
    15951564{
    1596 //printf("p_ExtGcd1:\n");
    1597 //p_Write(p, r); p_Write(q, r);
    15981565  assume((p != NULL) || (q != NULL)); 
    15991566  poly a = p; poly b = q; BOOLEAN aCorrespondsToP = TRUE;
  • libpolys/tests/polys_test.h

    r09ed77 rcfb500  
    209209{
    210210private:
     211  /* replaces n by n + c * var^exp (in-place);
     212     expects exp >= 0 */
     213  void plusTerm(number &n, int c, int exp, const coeffs cf)
     214  {
     215    number x = n_Par(1, cf);
     216    number pow = n_Init(1, cf);
     217    number temp;
     218    for (int i = 1; i <= exp; i++)
     219    {
     220      temp = n_Mult(pow, x, cf);
     221      n_Delete(&pow, cf);
     222      pow = temp;
     223    }
     224    n_Delete(&x, cf);
     225    number m = n_Init(c, cf);
     226    number t = n_Mult(m, pow, cf);
     227    n_Delete(&m, cf); n_Delete(&pow, cf);
     228    temp = n_Add(n, t, cf);
     229    n_Delete(&n, cf); n_Delete(&t, cf);
     230    n = temp;
     231  }
    211232  void checkInverse(number n, const coeffs cf)
    212233  {
     
    216237    number n2 = n_Mult(n, n1, cf);
    217238    clog << "(check: n * n^(-1) = "; n_Write(n2, cf);
     239    TS_ASSERT( n_IsOne(n2, cf) );
    218240    n_Delete(&n1, cf); n_Delete(&n2, cf);
    219241  }
     
    603625   
    604626    clog << endl
    605          << "Now let's compute some inverses in Q[a]/<a2+1>..."
     627         << "Now let's compute some inverses in Q[a]/<a^2+1>..."
     628         << endl;
     629   
     630    number u;
     631    u = NULL; plusTerm(u, 1, 1, cf); plusTerm(u, 1, 0, cf);  // a + 1
     632    checkInverse(u, cf); n_Delete(&u, cf);
     633    u = NULL; plusTerm(u, 1, 1, cf); plusTerm(u, -1, 0, cf); // a - 1
     634    checkInverse(u, cf); n_Delete(&u, cf);
     635    u = NULL; plusTerm(u, 1, 1, cf); plusTerm(u, 5, 0, cf);  // a + 5
     636    checkInverse(u, cf); n_Delete(&u, cf);
     637    u = NULL; plusTerm(u, 1, 1, cf); plusTerm(u, -5, 0, cf); // a - 5
     638    checkInverse(u, cf); n_Delete(&u, cf);
     639    u = NULL; plusTerm(u, 17, 1, cf); plusTerm(u, 5, 0, cf); // 17a + 5
     640    checkInverse(u, cf); n_Delete(&u, cf);
     641
     642    rDelete(s); // kills 'cf' and 'r' as well
     643  }
     644  void test_Q_Ext_b()
     645  {
     646    clog << "Start by creating Q[b]..." << endl;
     647
     648    char* n[] = {"b"};
     649    ring r = rDefault( 0, 1, n);   // Q[b]
     650    TS_ASSERT_DIFFERS( r, NULLp );
     651
     652    PrintRing(r);
     653
     654    TS_ASSERT( rField_is_Domain(r) );
     655    TS_ASSERT( rField_is_Q(r) );
     656   
     657    TS_ASSERT( !rField_is_Zp(r) );
     658    TS_ASSERT( !rField_is_Zp(r, 11) );
     659
     660    TS_ASSERT_EQUALS( rVar(r), 1);
     661
     662    poly minPoly = p_ISet(1, r);                    // minPoly = 1
     663    p_SetExp(minPoly, 1, 7, r); p_Setm(minPoly, r); // minPoly = b^7
     664    minPoly = p_Add_q(minPoly, p_ISet(17, r), r);   // minPoly = b^7 + 17
     665    ideal minIdeal = idInit(1);                     // minIdeal = < 0 >
     666    minIdeal->m[0] = minPoly;                       // minIdeal = < b^7 + 17 >
     667
     668    n_coeffType type = nRegister(n_Ext, naInitChar);
     669    TS_ASSERT(type == n_Ext);
     670
     671    ExtInfo extParam;
     672    extParam.r = r;
     673    extParam.i = minIdeal;
     674   
     675    clog << "Next create the extension field Q[b]/<b^7+17>..." << endl;
     676
     677    const coeffs cf = nInitChar(type, &extParam);   // Q[b]/<b^7+17>
     678   
     679    if( cf == NULL )
     680      TS_FAIL("Could not get needed coeff. domain");
     681
     682    TS_ASSERT_DIFFERS( cf->cfCoeffWrite, NULLp );
     683 
     684    if( cf->cfCoeffWrite != NULL )
     685    {
     686      clog << "Coeff-domain: "  << endl;
     687      n_CoeffWrite(cf); PrintLn();
     688    }
     689   
     690    // some tests for the coefficient field represented by cf:
     691    TestArithCf(cf);
     692    TestSumCf(cf, 10);
     693    TestSumCf(cf, 100);
     694    TestSumCf(cf, 101);
     695    TestSumCf(cf, 1001);
     696    TestSumCf(cf, 9000);
     697
     698    clog << "Finally create the polynomial ring (Q[b]/<b^7+17>)[u, v, w]..."
     699         << endl;
     700   
     701    char* m[] = {"u", "v", "w"};
     702    ring s = rDefault(cf, 3, m);   // (Q[b]/<b^7+17>)[u, v, w]
     703    TS_ASSERT_DIFFERS(s, NULLp);
     704
     705    PrintRing(s);
     706
     707    TS_ASSERT( rField_is_Domain(s) );
     708    TS_ASSERT( !rField_is_Q(s) );
     709    TS_ASSERT( !rField_is_Zp(s) );
     710    TS_ASSERT( !rField_is_Zp(s, 11) );
     711    TS_ASSERT( !rField_is_Zp(s, 13) );
     712    TS_ASSERT( !rField_is_GF(s) );
     713    TS_ASSERT( rField_is_Extension(s) );
     714    TS_ASSERT( !rField_is_GF(s, 25) );
     715    TS_ASSERT_EQUALS(rVar(s), 3);
     716
     717    Test(s);
     718   
     719    clog << endl
     720         << "Now let's compute some inverses in Q[b]/<b^7+17>..."
    606721         << endl;
    607722         
    608     number a = n_Par(1, cf);
    609     number n1 = n_Init(1, cf);
    610     number n5 = n_Init(5, cf);
    611     number n17 = n_Init(17, cf);
    612     number u; number v;
    613    
    614     u = n_Add(a, n1, cf); checkInverse(u, cf); n_Delete(&u, cf); // a + 1
    615     u = n_Sub(a, n1, cf); checkInverse(u, cf); n_Delete(&u, cf); // a - 1
    616     u = n_Add(a, n5, cf); checkInverse(u, cf); n_Delete(&u, cf); // a + 5
    617     u = n_Sub(a, n5, cf); checkInverse(u, cf); n_Delete(&u, cf); // a - 5
    618     v = n_Mult(a, n17, cf);
    619     u = n_Add(v, n5, cf); n_Delete(&v, cf);                      // 17a + 5
    620        checkInverse(u, cf); n_Delete(&u, cf);
    621    
    622     n_Delete(&a, cf); n_Delete(&n1, cf); n_Delete(&n5, cf);
    623     n_Delete(&n17, cf);
     723    number u;
     724    u = NULL; plusTerm(u, 1, 2, cf); plusTerm(u, 33, 0, cf);    // b^2 + 33
     725    checkInverse(u, cf); n_Delete(&u, cf);
     726    u = NULL; plusTerm(u, 1, 5, cf); plusTerm(u, -137, 0, cf);  // b^5 - 137
     727    checkInverse(u, cf); n_Delete(&u, cf);
     728
     729    clog << endl
     730         << "Now let's check a gcd computation in Q[b]..."
     731         << endl;
     732   
     733    number v;
     734    v = NULL; plusTerm(v, 1, 2, cf); plusTerm(v, 7, 1, cf);
     735    plusTerm(v, 1, 0, cf);        // b^2 + 7b + 1
     736    number w = n_Mult(v, v, cf);  // (b^2 + 7b + 1)^2
     737    number y = n_Mult(v, w, cf);  // (b^2 + 7b + 1)^3
     738    n_Delete(&v, cf);
     739    v = NULL; plusTerm(v, 2, 2, cf); plusTerm(v, -61, 1, cf);   // 2b^2 - 61b
     740    number z = n_Mult(w, v, cf);   // (b^2 + 7b + 1)^2 * (2b^2 - 61b)
     741    n_Delete(&v, cf);
     742   
     743    clog << "z = "; n_Write(z, cf);
     744    clog << "y = "; n_Write(y, cf);
     745    number theGcd = n_Gcd(z, y, cf);   // should yield w = (b^2 + 7b + 1)^2
     746    clog << "gcd(z, y) = "; n_Write(theGcd, cf);
     747   
     748    v = n_Sub(theGcd, w, cf);
     749    TS_ASSERT( v == NULL );
     750    n_Delete(&v, cf);
     751   
     752    clog << endl
     753         << "Now let's check an ext_gcd computation in Q[b]..."
     754         << endl;
     755         
     756    poly zFactor; poly yFactor;
     757    poly ppp = p_ExtGcd((poly)z, zFactor, (poly)y, yFactor, cf->algring);
     758    v = n_Sub(theGcd, (number)ppp, cf);
     759    TS_ASSERT( v == NULL );
     760    clog << "z = "; n_Write(z, cf);
     761    clog << "zFactor = "; p_Write(zFactor, cf->algring);
     762    clog << "y = "; n_Write(y, cf);
     763    clog << "yFactor = "; p_Write(yFactor, cf->algring);
     764    number v1 = n_Mult(z, (number)zFactor, cf);
     765    number v2 = n_Mult(y, (number)yFactor, cf);
     766    number v3 = n_Add(v1, v2, cf);
     767    clog << "z * zFactor + y * yFactor = "; n_Write(v3, cf);
     768    clog << "gcd(z, y) = "; p_Write(ppp, cf->algring);
     769    number v4 = n_Sub(v3, w, cf);
     770    TS_ASSERT( v4 == NULL );
     771   
     772    p_Delete(&ppp, cf->algring); p_Delete(&zFactor, cf->algring);
     773    p_Delete(&yFactor, cf->algring);
     774    n_Delete(&z, cf); n_Delete(&y, cf); n_Delete(&w, cf);
     775    n_Delete(&theGcd, cf); n_Delete(&v, cf); n_Delete(&v1, cf);
     776    n_Delete(&v2, cf); n_Delete(&v3, cf); n_Delete(&v4, cf);
    624777
    625778    rDelete(s); // kills 'cf' and 'r' as well
Note: See TracChangeset for help on using the changeset viewer.