воскресенье, 30 мая 2010 г.

Symbolic Integration

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

Классы функций одной переменной

Чем шире класс функции, тем сложнее брать интеграл. Например, интегрирование полиномов одной переменной
p(x) = an xn + ... + a1 x + a0
не составляет никакого труда. Уже некоторую сложность представляет интегрирование рациональной функции
r(x) = an xn + ... + a1 x + a0

bm xm + ... + b1 x + b0
.
Если степень числителя n меньше степени знаменателя m, и коэффициенты лежат в поле вещественных чисел, то интегрирование прозводится следующим образом: знаменатель раскладывается на множители вида (x - xi)ki и (x2 + bj x + cj)lj, где каждый xi — корень кратности ki, а bj и cj задают пару сопряжённых комплексных корней кратности lj. Затем методом неопределённых коэффициентов дробь r(x) раскладывается на сумму конечного количества дробей вида
Ai

(x-xi)ki
и Bj x+Cj

(x2 +bjx + cj)lj
Для каждого из трёх видов слагаемых (трёх, т.к. дробь с Bjx + Cj разбивается на отдельные слагаемые с Bjx и Cj) интеграл находится с помощью формул из справочника. В этих формулах встречаются логарифмы, радикалы, модули и арктангенсы, поэтому интеграл рациональной функции в общем случае не является рациональной функцией.

Если степень числителя n больше или равна степени знаменателя m, то мы можем разделить числитель на знаменатель и получить сумму полинома и «правильной» дроби. Если нас интересуют функции с комплексными коэффициентами, то всё делается точно так же, как и с вещественными, с тем отличием, что квадратичные члены из разложения знаменателя будут упразднены.

Кажется, что интегрирование рациональных функций — уже чисто техническая задача? Это не совсем так. Если степень многочлена больше 4, то его корни в общем случае не записываются в радикалах. Это доказал Абель в 1826 г. Например, вещественный корень уравнения «x5 + x + a = 0» невозможно записать, пользуясь арифметическими операциями и операцией извлечения корня.1

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

Алгебраическая функция y от переменной x задаётся неявным образом, через полином от x и y:
pn(x) yn + ... + p1(x) y + p0(x) = 0
Для примера преобразуем указанное выше уравнение 5-й степени в алгебраическую функцию y от переменной x:
y5+y+x=0
Алгебраические числа и функции тоже не являют собой предел выразительности математической записи. Хорошо известны трансцендентные числа (e, π) и трансцендентные функции (экспонента, логарифм, тригонометрические функции), которые невозможно представить через полиномиальные уравнения.

К классу элементарных функций относятся алгебраические функции, экспонента, логарифм и любые их конечные композиции. Подразумевается оперирование над полем комплексных чисел, что означает, что элементарны и все тригонометрические функции, т.к. они выражаются через комбинации eix и e-ix. Обратные тригонометрические функции, соответственно, выражаются через логарифмы. Показательная функция ax выражается через экспоненту и логарифм: ax = ex ln(a).

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

Интегрирование в теории

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

С алгебраическими функциями дело обстоит сложнее.

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

Это даёт следующую общую запись для интегрируемых алгебраических функций:
f(x) = v0'(x) + m
Σ
i=1
ci vi'(x)

vi(x)
(1)

Кроме того, Лаплас утверждал, что интеграл функции можно записать так, чтобы он содержал лишь те алгебраические величины, которые есть в исходном подынтегральном выражении.

Гипотезу Лапласа доказал Лиувилль в 1833 г., а ещё через два года он опубликовал доказательство аналогичного факта для некоторого обобщения класса алгебраических функций, т.н. элементарных функций по Лиувиллю. Они имеют форму алгебраических функций, но в роли переменных в них не только x, но ещё конечный набор логарифмов и экспонент от алгебраических функий от x:

y(x) = f(x, z1(x), ... zr(x)) — алгебраическая функция от x, z1(x), ... zr(x)
zi(x) = ln(qi(x)) или exp(qi(x)), где qi — алгебраическая функция от x.

Лиувилль показал, что если элементарная (по Лиувиллю) функция имеет элементарную первообразную, то эта первообразная представляется следующим образом:

f(x, z1(x), ... zr(x)) dx = y0(x, z1(x), ... zr(x)) + m
Σ
i=1
Ai ln(yi(x, z1(x), ... zr(x)))
Выражаясь менее формальным языком, первообразная функции y — это алгебраическая функция от того, что уже имеется в y, плюс конечное количество логарифмов (с коэффициентами) алгебраических функций того, что уже имеется в y. То есть, принцип Лапласа остаётся в силе.

Результаты Лиувилля, при всей своей значимости, совершенно не помогают в задаче нахождения первообразной, ибо процедуры приведения функции к интегрируемому виду автором предложено не было. Прошло более 80 лет до того, как английский математик Г. Х. Харди... предположил, что такой процедуры не может существовать, и ошибся.

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

Итак, дифференциальная теория Галуа. Ключевое понятие в ней — это т.н. элементарное расширение поля функций2, которое получается последовательным добавлением в него логарифмов, экспонент и алгебраических функций:
  • ln(f(x)), где f ∈ K
  • exp(f(x)), где f ∈ K
  • алгебраической функции, полученной уравнением вида P(f1(x), ... fn(x))=0, где f1, ... fn ∈ K.
(за K считается поле, расширенное предыдущим добавлением). Нетрудно видеть, что любая элементарная функция принадлежит какому-либо элементарному расширению поля рациональных функций.

Ключевая теорема дифференциальной теории Галуа является обобщением теоремы Лиувилля и говорит о том, что если функция f из поля F имеет первообразную g в поле G — элементарном расширении F, то f записывается в форме (1), где vi ∈ F, а g соответственно содержит функцию из F и логарифмы функций из F с весами.

Если строить поле F «специально» под функцию f, т.е. подбирать минимальное расширение поля рациональных функций, содержащее F, то станет очевидно, что принцип Лапласа снова не нарушается: первообразная функции f содержит то, что уже имеется в f, а также логарифмы того, что уже имеется в f.

Роберт Риш из Калифорнийского университета в 1969-1970 годах (Харди к тому времени более двадцати лет как не было в живых) опубликовал алгоритм, приводящий любую элементарную функцию к виду (1), или определяющий, что такое приведение невозможно (и следовательно, интеграла как элементарной функции нет). На этом вопрос отнюдь не был закрыт; напротив, началось самое интересное.

Интегрирование на практике

Алгоритм Риша, вопреки распространённому мнению, являлся на время своего опубликования алгоритмом интегрирования лишь чисто трансцендентных функций. Для алгебраических составляющих Риш привёл в общих чертах доказательство существования алгоритма, но чтобы сделать из него алгоритм, потребовался ещё не один десяток лет.

Первая реализация алгоритма Риша была выполнена Джоэлом Мозесом в рамках знаменитого «Project MAC» в MIT в 1971 г. Программа под названием SIN интегрировала чисто трансцендентные функции.

Дэвенпорт в 1981 г., основываясь на работе Риша и некоторых глубоких результатах дифференциальной алгебры и комплексного анализа, разработал алгоритм интегрирования чисто алгебраических функций и реализовал его в известной среде символьных вычислений REDUCE-2.

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

Алгоритм Дэвенпорта ещё и обладал большой вычислительной сложностью, что его автор признавал и выражал надежду, что алгоритм можно упростить. Так и произошло: Барри Трагер из MIT в 1984 г. внёс серьёзные улучшения в алгоритм Дэвенпорта. Обновлённый алгоритм обладал гораздо большей практической ценностью и был реализован в математических программах Axiom и Maple.

Решающий шаг к практическому решению вопроса сделал Мануэль Бронштейн в 1990 г, обобщив алгоритм Трагера на произвольные элементарные функции.

В 1998 г. Бронштейн написал монографию по символьному интегрированию, в которой понятным языком изложил лучшие достижения в этой области, начиная с теоремы Лиувилля и заканчивая собственными результатами.

Бронштейн был одним из ведущих разработчиков Axiom. Вместе с Трагером он реализовал в Axiom свой алгоритм, но лишь частично. В 2005 году Бронштейн скончался, спустя два года после смерти основателя Axiom Ричарда Дженкса.

При своей неполноте, имплементация Бронштейна — наиболее полная из имеющихся в программах символьных вычислений. Автор отмечал, что при возникновении нереализованного случая Axiom выдаёт ошибку, в отличие от процедур интегрирования Maple и Mathematica. То есть, если Axiom обнаруживает, что интеграл от данной функции не является элементарной функцией, то это правда. Если Axiom выдаёт ошибку, то интеграл может быть или не быть элементарной функцией.

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


x4 + 10 x2 - 96 x - 71
равна
- 1

8
ln((x6+15 x4-80 x3+27 x2-528 x+781)

x4+10 x2-96 x-71
- x8 - 20 x6 + 128 x5 - 54 x4 + 1408 x3 - 3124 x2 - 10001),
однако в таком виде её находит только Axiom. Прочие пакеты — Mathematica, Maple, Matlab, Maxima, Reduce — в лучшем случае выдают громоздкий ответ с участием эллиптических интегралов, т.е. неспособны выдать ответ в виде элементарной функции.3

Первообразную функции earcsin(x), равную ½(√ 1-x2 +x)earcsin(x), способны найти Mathematica, Axiom и Reduce.

Неопределённый интеграл ∫√arctan(x) dx невозможно найти (или определить его отсутствие) ни в одном из перечисленных пакетов; при этом Axiom выводит наиболее осмысленное сообщение: implementation incomplete (constant residues).

Нерешённые вопросы

Полная имплементация алгоритма Риша-Трагера-Бронштейна ждёт своего автора. Однако, и это ещё не всё: в теории символьного интегрирования тоже остались открытые вопросы.

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

Дело в том, что в процессе работы алгоритм Риша выстраивает «башню» элементарных расширений поля рациональных функций, с тем чтобы «верхушка» башни содержала подынтегральное выражение. При этом не всегда понятно, готова ли уже «башня» или нужно добавить в неё ещё один логарифм, экспоненту, или алгебраическую функцию. Например, если в подынтегральном выражении встретилась экспонента ex, она добавляется в поле (происходит элементарное расширение поля функцией ex). Если при дальнейшем анализе в подынтегральном выражении встречается ex2, то эту функцию также следует добавить в поле, вырастив таким образом «башню» до двух этажей. Но если после этого попадётся функция ex+x2, то алгоритм должен определить, что она принадлежит сформированному на данный момент элементарному расширению, и не наращивать «башню». Да, в описанном случае это несложно сделать; чуть сложнее определить, например, на языке элементарных функций, что arctan(x) = -arctan(1/x). Но в общем случае эта задача очень сложна4, и нынешние реализации алгоритма Риша используют лишь эвристики для её решения, которые по определению ненадёжны.

Как указывает Дэвенпорт и следом за ним Бронштейн, при реализации алгоритма можно «сделать вид», что нам известен т.н. базис трансцендентности для поля констант, то есть например считать e и √ 2 алгебраически независимыми. Это значит, что мы принимаем на веру, что нет никакого полинома от e и √ 2 , который был бы равен 0. Такое решение нужно принимать о каждой новой константе, появляющейся в ходе работы алгоритма, которую эвристика не позволяет связать алгебраическим выражением с ранее учтёнными константами. Если полученный базис трансцендентности окажется неверным, это может отразиться на корректности выданного ответа в том случае, если ответ гласит, что «интеграл не является элементарной функцией» (да, это касается и Axiom, о которой выше было сказано противоположное). Если же интеграл найден, его можно считать правильным.

Известным частным случаем задачи об алгебраической независимости является т.н. Constant Problem: определение того, тождественно ли нулю некоторое математическое выражение. Задача тривиально решается для многочленов, но для более широкого класса функций, как показал Даниэль Ричардсон в 1968 г., эта задача неразрешима.

Не всё, однако, так безнадёжно: Ричардсон в набор «элементарных расширений» включил, помимо экспоненты и логарифма, функцию модуля. В оригинале теорема Ричардсона звучит так:

Если есть функция f, получающаяся с помощью сложения, умножения и композиции из синуса, экспоненты и модуля, а коэффициенты f используют рациональные числа, π и ln(2), то задача проверки утверждения «f тождественно равна нулю» алгоритмически неразрешима.

Функция модуля существенна в доказательстве Ричардсона, и есть основания полагать (да и сам автор теоремы надеется), что без неё задача Constant Problem всё таки разрешима, что даёт надежду на создание абсолютно честного алгоритма символьного интегрирования.

И последнее: класс элементарных функций кажется несколько надуманным, не так ли? Например, функция ошибки
erf(x) = 2



π

x

0
e-t2 dt
не является элементарной, как доказал Лиувилль. Но чем она хуже экспоненты и логарифма, которые мы относим к элементарным функциям? На компьютере она может быть вычислена с произвольной точностью, и даже, как экспонента и логарифм, обладает элементарной производной. Эту функцию вполне можно включить в «надэлементарные» расширения функций, так же как и эллиптические интегралы, ∫xxdx и многие другие. Тема интегрируемости в «надэлементарных» функциях была затронута ещё автором SIN в 1970-х годах, и ныне представляет обширное поле для исследований.

Благодарности

Автор признателен А. В. Демьянову, П. А. Мозоляко и Н. Н. Васильеву за вычитку черновиков этой статьи и высказанные замечания.

Литература

J. Liouville. Premier memoire sur la determination des integrales dont la valeur est algebrique. Journal de l’Ecole Polytechnique, 14:124–148, 1833

J. Liouville. Second memoire sur la détermination des integrales dont la valeur est algebrique. Journal de l’Ecole Polytechnique, 14:149–193, 1833

J. Liouville. Mèmoire sur l’intégration d’une classe de fonctions transcendantes. J. Reine Angew. Math. Bd. 13, p. 93-118. (1835)

G. H. Hardy. The Integration of Functions of a Single Variable, 1916.

R. H. Risch (1970). The solution of the Problem of Integration in Finite Terms. Bulletin of the American Mathematical Society 76 (3): 605–608.

Joel Moses. Symbolic integration the stormy decade. Proceedings of the second ACM symposium on Symbolic and algebraic manipulation (1971) 427-440.

Дж. Дэвенпорт. Интегрирование алгебраических функций, Москва «Мир» 1985.

B. Trager. Integration of algebraic functions. PhD Thesis, MIT, 1984

M. Bronstein. Integration of elementary functions. Journal of Symbolic Computation, Volume 9, Issue 2 (1990). 117-173

M. Bronstein. Symbolic Integration Tutorial. ISSAC'98, Rostock.

Richardson, D. Some Unsolvable Problems Involving Elementary Functions of a Real Variable. J. Symbolic Logic 33, 514-520, 1968.

Сноски:

1Кстати, это число называется корнем Бринга Br(a), и с его помощью можно записать все комплексные корни более общего уравнения «x5 + px + q = 0».
2Формально говоря, дифференциальная теория Галуа оперирует дифференциальными полями характеристики 0, но мы в эти термины погружаться не станем и будем считать по-прежнему, что имеем дело с различными классами функций.
3Проверялись следующие версии программ: Axiom 2010, Mathematica 7, Maple 13, Matlab R2007, Maxima 5.21, Reduce 2009.
4Например, до сих пор никто не знает, является ли e + π рациональным числом (хотя известно, что e и π по отдельности — трансцентентные числа). Известно только равенство ei π = -1 и бесконечная дробь Рамануджана, связывающая e и π.