Вернуться к списку примеров

Многочлены с рациональными коэффициентами

по материалам курсовой работы О. Пиценко


Хорошо известно, что если рациональное число p/q является корнем уравнения с целыми коэффициентами anxn+ an-1xn-1+ ... a1x+ a0=0, то a0 делится на p и an делится на q. Это свойство дает возможность по коэффициентам a0 и an данного уравнения найти все рациональные числа, которые могут быть корнями этого уравнения. Подставляя затем все эти числа в уравнение, можно установить, какие из них действительно являются его корнями. При этом переборе обычно используют еще одно необходимое условие: для того, чтобы рациональное число p/q являлось корнем уравнения f(x)=0 с целыми коэффициентами, необходимо, чтобы значение функции f(1) делилось на p-q, а значение f(-1) делилось на p+q.

Имея уравнение с рациональными коэффициентами, легко получить эквивалентное ему уравнение с целыми коэффициентами. Следующая функция использует перечисленные выше факты для вычисления рациональных корней многочлена с рациональными коэффициентами, при этом выводя на печать необходимые промежуточные вычисления.

RationalRootsOfPolynomial:=function(f)
local x, a0, diva0, an, g, divan, roots, p, q, x0, val, a, i, c, k,
      fplus1, fminus1, t, d;
Print("f(x) = ", f,"\n");
x:=IndeterminateOfUnivariateRationalFunction(f);
# Создаем пустой список для хранения в нем корней многочлена
roots:=[]; 
# Oпределяем коэффициенты 
a:=PolynomialCoefficientsOfPolynomial(f,x);
# Переводим коэффициенты в числа 
a:=List(a, i -> Value(i,0) ); 
a:=List(a, i -> DenominatorRat(i) );
# Находим общий делитель знаменателей всех коэффициентов
a:=Lcm(a); 
# Домножаем многочлен на общий знаменатель
f:=f*a; 
if a>1 then
  Print("f(x) = 1/", a, " * (", f, ")\n");
fi;
a0:=Value(f,0); 
# Находим свободный член a0 и проверяем его равенство нулю. 
# Если a0=0, то делим многочлен на x и добавляем к списку roots
# корень уравнения 0 столько раз, какова его кратность
if a0=0 then
  k:=PositionProperty( PolynomialCoefficientsOfPolynomial(f,x), 
                       c -> not IsZero(c) ) - 1;
  f:=Quotient( f, x^k );
  AddSet(roots, 0);
  a0:=Value(f,0);
  if a=1 then
    Print("f(x) = x^", k, " * (", f,")\n");
  else
    Print("f(x) = 1/", a, " * x^", k, " * (", f,")\n");
  fi;
fi;
Print("a_0 = ", a0, "\n");
# Находим делители a0
diva0:=DivisorsInt(a0); 
Print("diva0 = ", diva0, "\n");
an:=LeadingCoefficient(f);
Print("a_n = ", an, "\n");
# Находим делители an
divan:=DivisorsInt(an); 
Print("divan = ", divan, "\n");
# Находим значение многочлена в точкe 1
fplus1:=Value(f,1); 
# Находим значение многочлена в точкe -1
fminus1:=Value(f,-1); 
Print(" p      q    (p-q)|f(1)      (p+q)|f(-1)    f(p/q)", "\n");
# Проверяем делимость f(1) на (p-q) и f(-1) на (p+q). В случае положительного 
# результата подстановкой проверяем, является ли число p/q корнем данного многочлена.
# Если число p/q является корнем многочлена, то добавляем его в список roots
for t in diva0 do
  for q in divan do
    for p in [t, -t] do
      x0:=p/q;
      Print(" ",p,"      ",q, "\c");
      d:=p-q;
      if d<>0 then
        if (fplus1 mod (p-q)) = 0 then
          Print("       +       ", "\c");
        else
          Print("       -       ", "\n");
          continue;
        fi;
      else
        Print( " ", "\c");
      fi;
      d:=p+q;
      if p+q<>0 then
        if (fminus1 mod (p+q)) = 0 then
          Print("       +       ", "\c");
        else
          Print("       -       ", "\n");
          continue;
        fi;
      else
        Print( " ", "\c");
      fi;
      val:=Value(f,x0);
      Print(" f(", x0, ") = ", val, "\n");
      if val=0 then
        AddSet(roots, x0);
        # Найдя корень многочлена x0, делим исходный многочлен на (x-x0)
        f:=Quotient(f,x-x0); 
      fi;
    od;
  od;
od;
return roots;
end;

Приведем пример работы с данной функцией (в нашу задачу не входила дальнейшая автоматизация представления информации, поэтому доработку программы для выравнивания данных мы оставляем читателю в качестве упражнения).


Пример 1.

gap> Q:=Rationals;
Rationals
gap> x:=Indeterminate(Q,"x");
x
gap> f:=(x^2-4/9)*x^2;
x^4-4/9*x^2
gap> RationalRootsOfPolynomial(f);
f(x) = x^4-4/9*x^2
f(x) = 1/9 * (9*x^4-4*x^2)
f(x) = 1/9 * x^2 * (9*x^2-4)
a_0 = -4
diva0 = [ 1, 2, 4 ]
a_n = 9
divan = [ 1, 3, 9 ]
 p      q    (p-q)|f(1)      (p+q)|f(-1)    f(p/q)
 1      1        -       
 -1      1       -       
 1      3       -       
 -1      3       -       
 1      9       -       
 -1      9       -       
 2      1       +              -       
 -2      1       -       
 2      3       +              +        f(2/3) = 0
 -2      3       +              +        f(-2/3) = 0
 2      9       -       
 -2      9       -       
 4      1       -       
 -4      1       +              -       
 4      3       +              -       
 -4      3       -       
 4      9       +              -       
 -4      9       -       
[ -2/3, 0, 2/3 ]
gap>
Пример 2.

gap> f:=(x^2-1/4)*(x+5/3);
x^3+5/3*x^2-1/4*x-5/12
gap> RationalRootsOfPolynomial(f);
f(x) = x^3+5/3*x^2-1/4*x-5/12
f(x) = 1/12 * (12*x^3+20*x^2-3*x-5)
a_0 = -5
diva0 = [ 1, 5 ]
a_n = 12
divan = [ 1, 2, 3, 4, 6, 12 ]
 p      q    (p-q)|f(1)      (p+q)|f(-1)    f(p/q)
 1      1        +        f(1) = 24
 -1      1       +         f(-1) = 6
 1      2       +              +        f(1/2) = 0
 -1      2       +              +        f(-1/2) = 0
 1      3       +              -       
 -1      3       +              +        f(-1/3) = 16
 1      4       +              -       
 -1      4       -       
 1      6       -       
 -1      6       -       
 1      12       -       
 -1      12       -       
 5      1       +              +        f(5) = 80
 -5      1       +              -       
 5      2       +              -       
 -5      2       -       
 5      3       +              -       
 -5      3       +              +        f(-5/3) = 0
 5      4       +              -       
 -5      4       -       
 5      6       +              -       
 -5      6       -       
 5      12       -       
 -5      12       -       
[ -5/3, -1/2, 1/2 ]
gap>
Созданную нами функцию RationalRootsOfPolynomial удобно использовать для проверки заданий на нахождение рациональных корней многочленов в курсе алгебры и теории чисел, однако она не является универсальной. В частности, эта функция не учитывает кратность корней. Поэтому мы разработаем еще одну функцию MyRootsOfUPol, которая будет возвращать рациональные корни многочленов с учетом их кратности:

MyRootsOfUPol:=function(f) 
local x, a0, diva0, an, g, divan, roots, p, q, x0, val, a, i, c, k, 
      fplus1, fminus1, t, d, var; 
if IsZero(f) then  
  return [0]; 
elif DegreeOfUnivariateLaurentPolynomial(f)=0 then 
  return []; 
fi; 
x:=IndeterminateOfUnivariateRationalFunction(f); 
roots:=[]; 
a:=PolynomialCoefficientsOfPolynomial(f,x); 
a:=List(a, i -> Value(i,0) ); 
a:=List(a, i -> DenominatorRat(i) ); 
a:=Lcm(a);  
if a<>1 then 
  f:=f*a; 
fi;  
a0:=Value(f,0); 
if a0=0 then 
  k:=PositionProperty( PolynomialCoefficientsOfPolynomial(f,x), 
                       c -> not IsZero(c) ) - 1; 
  f:=Quotient( f, x^k ); 
  Append(roots, List([1..k], i -> 0) ); 
  a0:=Value(f,0); 
fi; 
diva0:=DivisorsInt(a0); 
an:=LeadingCoefficient(f); 
divan:=DivisorsInt(an); 
fplus1:=Value(f,1); 
fminus1:=Value(f,-1); 
var:=[]; 
for t in diva0 do 
  for q in divan do 
    for p in [t, -t] do 
      AddSet(var, p/q); 
    od; 
  od; 
od; 
for x0 in var do 
  p:=NumeratorRat(x0); 
  q:=DenominatorRat(x0); 
  d:=p-q; 
  if d<>0 then  
    if (fplus1 mod (p-q)) <> 0 then 
      continue; 
    fi; 
  fi; 
  d:=p+q; 
  if p+q<>0 then  
    if (fminus1 mod (p+q)) <> 0 then 
      continue; 
    fi; 
  fi; 
  val:=Value(f,x0); 
  if val=0 then  
    Add(roots, x0); 
    if DegreeOfUnivariateLaurentPolynomial(f)>0 then 
      f:=Quotient(f,x-x0); 
      Append(roots, MyRootsOfUPol(f)); 
      return roots; 
    fi;  
  fi;  
od; 
return roots; 
end; 

Алгоритм данной функции практически совпадает с предыдущей, однако она работает быстрее, так как не тратит время на вывод промежуточных вычислений. Интересно сравнить производительность функции MyRootsOfUPol со стандартной функцией системы GAP RootsOfUPol, а также проверить правильность расчетов сравнением результатов, возвращаемых обеими функциями. Зададим многочлен десятой степени, случайным образом сформировав его коэффициенты, и найдем его рациональные корни двумя способами:

gap> f:=Product(List([1..10],i->x-Random(Rationals)));
x^10-17/3*x^9+43/16*x^8+247/16*x^7+197/64*x^6-279/64*x^5-77/64*x^4+41/192*x^3+1/16*x^2
gap> RootsOfUPol(f);               
[ 4, 3, 1/2, 1/4, 0, 0, -1/4, -1/3, -1/2, -1 ]
gap> time;
283
gap> MyRootsOfUPol(f);
[ 0, 0, -1, -1/2, -1/3, -1/4, 1/4, 1/2, 3, 4 ]
gap> time;
9
gap> SortedList(MyRootsOfUPol(f))=SortedList(RootsOfUPol(f));
true
gap>
Из приведенного выше примера видно, что результаты расчетов совпадают. Как показывает этот и аналогичные примеры, при небольшой степени многочлена функция MyRootsOfUPol находит рациональные корни многочлена быстрее, чем стандартная функция RootsOfUPol. Возьмем теперь более сложный многочлен, например, многочлен 30-й степени и сравним результаты и время работы обеих функций:

gap> f:=Product(List([1..30],i->x-Random(Rationals)));
x^30+1051/60*x^29+141331/1200*x^28+3082613/8640*x^27+58522771/259200*x^26-436162163/259200*x^2\
5-35773349/7200*x^24-7211778329/2073600*x^23+33235831763/4147200*x^22+125514755/6912*x^21+7199\
839177/1036800*x^20-9780069169/518400*x^19-34659602443/1382400*x^18-3695473133/2073600*x^17+98\
3525023/51840*x^16+9485611669/691200*x^15-9502389071/4147200*x^14-7894762583/1036800*x^13-1052\
660579/345600*x^12+190257823/207360*x^11+4904751197/4147200*x^10+216549797/691200*x^9-1813583/\
28800*x^8-595997/9600*x^7-54551/3200*x^6-3513/1600*x^5-9/80*x^4
gap> RootsOfUPol(f);
[ 2, 1, 1, 1, 2/3, 2/3, 3/5, 0, 0, 0, 0, -1/5, -1/4, -1/2, -1/2, -1/2, -1/2, -3/4, -3/4, -1, 
  -1, -1, -1, -1, -4/3, -3/2, -5/3, -2, -3, -6 ]
gap> time;
102
gap> MyRootsOfUPol(f);
[ 0, 0, 0, 0, -6, -3, -2, -5/3, -3/2, -4/3, -1, -1, -1, -1, -1, -3/4, -3/4, -1/2, -1/2, -1/2, 
  -1/2, -1/4, -1/5, 3/5, 2/3, 2/3, 1, 1, 1, 2 ]
gap> time;
993
gap> SortedList(MyRootsOfUPol(f))=SortedList(RootsOfUPol(f));
true
gap> 
Данный пример показывает, что для сложных многочленов стандартная функция RootsOfUPol все-таки находит рациональные корни значительно быстрее. Таким образом, наша разработка не может претендовать на замещение стандартной функции системы. Однако, в отличие от RootsOfUPol, использующей поле разложения многочлена, функции, предлагаемые нами, используют минимальный набор теоретических сведений из стандартного курса алгебры и теории чисел, и могут являться полезным учебным примером, показывающим, что нужно на самом деле сделать для того, чтобы запрограммировать алгоритм, базирующийся всего лишь на на нескольких простых утверждениях.


Вернуться к списку примеров