Maxima

Материал из свободной русской энциклопедии «Традиция»
Перейти к навигации Перейти к поиску

Maxima — система компьютерной алгебры, реализованная на языке Lisp.

Maxima имеет широчайший набор средств для проведения аналитических вычислений, численных вычислений и построения графиков. По набору возможностей система близка к таким коммерческим системам как Maple и Mathematica. В то же время она обладает высочайшей степенью переносимости. Это единственная из существующих систем аналитических вычислений, которая может работать на всех современных операционных системах, на компьютерах, начиная от самых мощных, заканчивая наладонниками.

История[править | править код]

Maxima произошла от системы «Macsyma», разрабатывавшейся в MIT (Massachusetts Institute of Technology) с 1968 по 1982.

Вариант продукта (известный как DOE Macsyma) сопровождался профессором Уильямом Шелтером в Техасском Университете с 1982 до своей смерти в 2001.

В 1998, Шелтер получил от Департамента энергии разрешение опубликовать исходный код DOE Macsyma под лицензией GPL, и в 2000 он создал проект «Maxima» на SourceForge для поддержания и развития «DOE Macsyma», переименованного в Maxima.

Введение[править | править код]

Обзор[править | править код]

Maxima может использоваться как мощный калькулятор:

(%i13) 144*17 - 9;
(%o13) 2439


Maxima поддерживает точные вычисления с большими числами, размер которых ограничен только ресурсами вашего компьютера.

Например, вычислим <m>144^{25}</m>, а затем возьмем корень 25-степени из результата:

(%i14) 144^25;
(%o14) 910043815000214977332758527534256632492715260325658624
(%i15) %^(1/25);
(%o15) 144

Знак «%» ссылается на специльную переменную, хранящую результат последней операции. Знак «^» — оператор возведения в степень.

Однако основное назначение компьютерной алгебры — не численные вычисления, а скорее, символьные операции. Поиграем с полиномом <m>(x+2y)^4</m>:

expand((x + 2*y)^4); → <m>16\,y^4+32\,x\,y^3+24\,x^2\,y^2+8\,x^3\,y+x^4</m>
factor(%);           → <m>\left(2\,y+x\right)^4</m>


Вычисление производных:

diff(sin(x)*cos(x), x);     → <m>\cos ^2x-\sin ^2x</m>
trigsimp(%);                → <m>2\,\cos ^2x-1</m>
diff(%, x);                 → <m>-4\,\cos x\,\sin x</m>
diff( sin(x)*cos(x), x, 2); → <m>-4\,\cos x\,\sin x</m>

и первообразных (неопределенных интегралов):

integrate((x + 1)/(x^3 - 8), x); → <m>-{{\log \left(x^2+2\,x+4\right)}\over{8}}+{{\arctan \left({{2\,x+2
}\over{2\,\sqrt{3}}}\right)}\over{4\,\sqrt{3}}}+{{\log \left(x-2
\right)}\over{4}}</m>

Maxima может преобразовывать тригонометрические выражения в канонические конечные тригонометрические ряды (Fourier sums):

trigreduce (sin(x)^5); → <m>{{\sin \left(5\,x\right)-5\,\sin \left(3\,x\right)+10\,\sin x
}\over{16}}</m>


Первые шаги[править | править код]

Maxima ориентирована на минималистический интерфейс командной строки — т. е. сеанс работы состоит из последовательных команд пользователя и полученных на них ответов.

Команда — это некое выражение, возможно многострочное, оканчивающееся точкой с запятой «;». Например на команду:

(%i1) 14*(7 + 3);

будет получен ответ:

(%o1)  140

[[#ref_{{{1}}}|]] Если вы забудете ввести «;» — Maxima будет терпеливо (т. е. без каких-либо сообщений об ошибках) продолжать ввод ввод, пока вы ее не введете «;».

Шаблон:Tip Можно вместо точки с запятой завершить команду знаком доллара «$», тогда вывод этой команды будет заглушён, то есть, команда будет обработана и результат будет вычислен, но не выведен на экран.

Введем какое-нибудь выражение с переменной:

(%i2) (x + 2)*(x - 3);

Ответ Maxima содержит то же выражение, в неупрощенной форме:

(%o2)  (x - 3) (x + 2)

[[#ref_{{{1}}}|]]  Maxima выдает выражение-результаты, используя минималистические средства — текстовую отрисовку формул. Так в ответе «(%o2)» мы видим, что за счет исключения «*» из ответа, выражение выглядит «математическим», а не программистким.

Теперь раскроем скобки:

(%i3) expand(%);
				   2
(%o3) 				  x  - x - 6

Т.е. человек вполне способен прочесть эту ASCII-формулу, однако, есть возможность получить эту формулу (c помощью функции «tex») в формате TeX:

(%i38) tex( expand( (x + 2)*(x - 3) ) );
$x^2-x-6$
(%o38) false 

Что, если отрисовать через TeX будет выглядеть так: <m>x^2-x-6</m>.

[[#ref_{{{1}}}|]]  Пользуясь этой возможностью, далее в этом введении, все математические формулы, получаемые от Maxima мы будем отрисовывать через TeX, что будет гораздо удобней для чтения, чем ASCII-графика, особенно в случае «многоэтажных» формул.

Итак, возведение в степень:

x^3^4;        → <m>x^{81}</m>

Обратите внимание на приоритет операций: x^3^4 есть <m>x^{3^4}=x^{81} \ne (x^3)^4 = x^{3\cdot 4}= x^{12} </m>

Дифференцирование и интегрирование[править | править код]

Возьмем неопределенный интеграл для <m>{{x}\over{x^3+1}}</m>:

integrate(x/(x^3 + 1), x);  → <m>{{\log \left(x^2-x+1\right)}\over{6}}+{{\arctan \left({{2\,x-1
}\over{\sqrt{3}}}\right)}\over{\sqrt{3}}}-{{\log \left(x+1\right)
}\over{3}}</m>

Продифференцируем результат для проверки:

diff(%, x); → <m>{{2}\over{3\,\left({{\left(2\,x-1\right)^2}\over{3}}+1\right)}}+{{2
\,x-1}\over{6\,\left(x^2-x+1\right)}}-{{1}\over{3\,\left(x+1\right)
}}</m>

И чтобы, не пугаться, упростим сумму рациональных выражений:

ratsimp(%); → <m>{{x}\over{x^3+1}}</m>

Повторим тоже самое для <m>e^{a\,x}\,\cos x\,\sin x</m>:

integrate(exp(a*x)*sin(x)*cos(x), x); → <m>{{e^{a\,x}\,\left(a\,\sin \left(2\,x\right)-2\,\cos \left(2\,x
\right)\right)}\over{2\,\left(a^2+4\right)}}</m>
 diff(%, x); → <m>{{a\,e^{a\,x}\,\left(a\,\sin \left(2\,x\right)-2\,\cos \left(2\,x
\right)\right)}\over{2\,\left(a^2+4\right)}}+{{e^{a\,x}\,\left(4\,
\sin \left(2\,x\right)+2\,a\,\cos \left(2\,x\right)\right)}\over{2\,
\left(a^2+4\right)}}</m>
 ratsimp(%); → <m>{{e^{a\,x}\,\sin \left(2\,x\right)}\over{2}}</m>
trigexpand(%); → <m>e^{a\,x}\,\cos x\,\sin x</m>

Если Maxima не смогла взять интеграл — она возвращает интеграл в «исходной форме»:

(%i4) integrate(1/((x-3)^4+1/2), x);
                        /
                        [      1
(%o4)                   I ------------ dx
                        ]        4   1
                        / (x — 3)  + -
                                     2

Шаблон:Tip Разумеется, в TeX-форматировании интеграл будет выглядеть красивее: <m>\int {{{1}\over{\left(x-3\right)^4+{{1}\over{2}}}}}{\;dx}</m>

В таком случае, можно помочь системе, выполнив замену переменных:

changevar (%, x - 3 - y ,y ,x); → <m>2\,\int {{{1}\over{2\,y^4+1}}}{\;dy}</m>

а затем проинтегрировав:

ev (%, integrate) 

<m>2\,\left({{\log \left(\sqrt{2}\,y^2+2^{{{3}\over{4}}}\,y+1\right)

}\over{4\,2^{{{3}\over{4}}}}}-{{\log \left(\sqrt{2}\,y^2-2^{{{3
}\over{4}}}\,y+1\right)}\over{4\,2^{{{3}\over{4}}}}}+{{\arctan 
\left({{2\,\sqrt{2}\,y+2^{{{3}\over{4}}}}\over{2^{{{3}\over{4}}}}}
\right)}\over{2\,2^{{{3}\over{4}}}}}+{{\arctan \left({{2\,\sqrt{2}\,
y-2^{{{3}\over{4}}}}\over{2^{{{3}\over{4}}}}}\right)}\over{2\,2^{{{3
}\over{4}}}}}\right)</m>

выполним обратную подстановку

sfx: %, y=x-3;

<m>2\,\left({{\log \left(\sqrt{2}\,\left(x-3\right)^2+2^{{{3}\over{4}}

}\,\left(x-3\right)+1\right)}\over{4\,2^{{{3}\over{4}}}}}-{{\log 
\left(\sqrt{2}\,\left(x-3\right)^2-2^{{{3}\over{4}}}\,\left(x-3
\right)+1\right)}\over{4\,2^{{{3}\over{4}}}}}+{{\arctan \left({{2\,
\sqrt{2}\,\left(x-3\right)+2^{{{3}\over{4}}}}\over{2^{{{3}\over{4}}}
}}\right)}\over{2\,2^{{{3}\over{4}}}}}+{{\arctan \left({{2\,\sqrt{2}
\,\left(x-3\right)-2^{{{3}\over{4}}}}\over{2^{{{3}\over{4}}}}}
\right)}\over{2\,2^{{{3}\over{4}}}}}\right)</m>

Т.к. подынтегральное выражение — всюду определенная положительная функция, мы можем вычислить определенный интеграл <m> 0 1 1 ( x 3 ) 4 + 1 2 d x \int_{0}^{1}{{{1}\over{\left(x-3\right)^4+{{1}\over{2}}}}\;dx} </m>, используя вычисленную первообразную «sfx»:

ratsimp(subst (1, x, sfx) - subst(0, x, sfx));

<m>{{2\,2^{{{1}\over{4}}}\,\arctan \left({{2^{{{3}\over{4}}}+6\,\sqrt{

2}}\over{2^{{{3}\over{4}}}}}\right)-2\,2^{{{1}\over{4}}}\,\arctan 
\left({{2^{{{3}\over{4}}}+4\,\sqrt{2}}\over{2^{{{3}\over{4}}}}}
\right)+2^{{{1}\over{4}}}\,\log \left(3\,2^{{{3}\over{4}}}+9\,\sqrt{
2}+1\right)-2^{{{1}\over{4}}}\,\log \left(2\,2^{{{3}\over{4}}}+4\,
\sqrt{2}+1\right)+2^{{{1}\over{4}}}\,\log \left(-2\,2^{{{3}\over{4}}
}+4\,\sqrt{2}+1\right)-2^{{{1}\over{4}}}\,\log \left(-3\,2^{{{3
}\over{4}}}+9\,\sqrt{2}+1\right)+2\,2^{{{1}\over{4}}}\,\arctan 
\left({{6\,\sqrt{2}-2^{{{3}\over{4}}}}\over{2^{{{3}\over{4}}}}}
\right)-2\,2^{{{1}\over{4}}}\,\arctan \left({{4\,\sqrt{2}-2^{{{3
}\over{4}}}}\over{2^{{{3}\over{4}}}}}\right)}\over{4}}</m>

или переходя к значениям в плавающей арифметике:

float(ratsimp(subst (1, x, sfx) - subst(0, x, sfx))); → 0.028806333852738

Для проверки выполним численное интегрирование методом Ромберга:

romberg(1/((x-3)^4+1/2), x,0, 1); → 0.028806333924554

Результаты аналитического и численного интегрирования совпадают с точностью 1e-10.

Проинтегрируем <m>\left(x^4+1\right)^{{{2}\over{3}}}\,\left(6\,x^5+7\,x^4-36\,x^3+18

\,x-21\right)</m>.
integrate(((6*x^5 + 7*x^4 - 36*x^3 + 18*x - 21)*(x^4 + 1)^(2/3)
       + (2*x^6 - 20*x^4 - 40*x^3 + 18*x^2 + 12)*(x^4 + 1)^(1/3))
        /(3*x^8 +6*x^4 + 3), x);

<m>{{\left(3\,x^2-7\,x+9\right)\,e^{{{2\,\log \left(x^4+1\right)

}\over{3}}}+\left(2\,x^3+4\,x+5\right)\,e^{{{\log \left(x^4+1\right)
}\over{3}}}}\over{x^4+1}}</m>

Осталось «сократить» экспоненты с логарифмами, перейдя к каноническому радикалу:

radcan(%);

<m>

{{\left(3\,x^2-7\,x+9\right)\,\left(x^4+1\right)^{{{2}\over{3}}}+
\left(2\,x^3+4\,x+5\right)\,\left(x^4+1\right)^{{{1}\over{3}}}
}\over{x^4+1}}

</m>

Решение полиномов[править | править код]

Нахождение нулей у полиномов от одной переменной — задача известная со школы. Поиграем с ней для демонстрации некоторых возможностей Maxima. Определим переменную poly:

poly: x^2 - x - 12; → <m>x^2-x-12</m>

Вычислим нули, используя функцию «solve»:

solutions: solve (poly=0, x); → [ x=-3, x=4 ]

Мы получили ответ в виде двух равенств (с левыми и правыми частями) — если их преобразовать в множители, то мы должны получить исходный полином, в виде их произведения.

Воспользуемся механизмом лямбда-функций (по сути, безымянных функций). Лямбда функция

lambda ([eq], lhs(eq) - rhs(eq))

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

map( lambda( [eq], lhs(eq) - rhs(eq)), solutions); → [x+3,x-4]

Перемножим этот список:

apply("*", %); → (x-4)(x+3)

И, раскрыв скобки, получаем исходный полином:

expand(%); → <m>x^2-x-12</m>

Тригонометрические преобразования[править | править код]

Чтобы преобразовать произведения синусов-косинусов, в триногометрический ряд Фурье, используется функция «trigreduce»:

trigreduce(sin(x)^3*cos(x)); → <m>{{2\,\sin \left(2\,x\right)-\sin \left(4\,x\right)}\over{8}}</m>

Теперь избавимся от произведений в аргументах тригонометрических функций:

trigexpand(%); → <m>{{4\,\cos x\,\sin ^3x-4\,\cos ^3x\,\sin x+4\,\cos x\,\sin x}\over{8}}</m>

И наконец, упрощаем выражение, используя функцию «trigsimp», знающую тождество <m>sin^2x + cos^2y = 1</m>:

trigsimp(%); → <m>\cos x\,\sin ^3x</m>

Упрощение тригонометрических выражений не всегда получается просто и быстро, как в приведенном примере, и в почтовых рассылках часто раздаются вопросы, «как мне упростить XXX?» или «почему Maxima не упрощает YYY?».

Рассмотрим пример:

result:integrate(9/cos(7*x)^2, x); → <m>{{18\,\sin \left(14\,x\right)}\over{7\,\sin ^2\left(14\,x\right)+7
\,\cos ^2\left(14\,x\right)+14\,\cos \left(14\,x\right)+7}}</m>

А пользуясь любым матсправочником с таблицей основных первообразных (например, справочником Бронштейна-Семендяева), легко получить более красивый результат: <m>\frac{9\tan(7x)}{7}</m>.

Попробуем упростить «влоб»:

simp1: trigsimp(result); → <m>{{9\,\sin \left(14\,x\right)}\over{7\,\cos \left(14\,x\right)+7}}</m>

Но применять дальше «trigexpand», для раскрытия sin(14*x) и cos(14*x) через sin(x) и cos(x) бесполезно — выражение сильно усложнится. Надо бы выразить sin(14*x) и cos(14*x) через sin(7*x) и cos(7*x), а для этого надо сделать подстановку, заменив 14*x на 2*y.

simp2: simp1, 14*x=2*y; → <m>{{9\,\sin \left(2\,y\right)}\over{7\,\cos \left(2\,y\right)+7}}</m>

И далее:

trigexpand(simp2); → <m>{{18\,\cos y\,\sin y}\over{7\,\left(\cos ^2y-\sin ^2y\right)+7}}</m>
trigsimp(%);       → <m>{{9\,\sin y}\over{7\,\cos y}}</m>

и наконец:

trigreduce(%);     → <m>{{9\,\tan y}\over{7}}</m>

возвращаемся к исходным переменным:

%, y = 7*x;        → <m>{{9\,\tan \left(7\,x\right)}\over{7}}</m>

Итак, мы показали на примере этой задачки, как важно использовать подстановки переменных при упрощении через «trigexpand». А вообще, это выражение можно было упростить одним действием:

trigrat(result);   → <m>{{9\,\sin \left(7\,x\right)}\over{7\,\cos \left(7\,x\right)}}</m>
trigreduce(%);     → <m>{{9\,\tan \left(7\,x\right)}\over{7}}</m> 

Проверим результат интегрирования:

deriv: diff(result, x);
trigrat(deriv);    → <m>{{18}\over{\cos \left(14\,x\right)+1}}</m>

Выглядит просто, но чтобы сравнить с исходным выражением, надо выразить всю «тригонометрию» через 7x.

%, 14*x = 2*y;   → <m>{{18}\over{\cos \left(2\,y\right)+1}}</m>
trigexpand(%);   → <m>{{18}\over{-\sin ^2y+\cos ^2y+1}}</m>
trigsimp(%);     → <m>{{9}\over{\cos ^2y}}</m>
%, y = 7*x;      → <m>{{9}\over{\cos ^2\left(7\,x\right)}}</m>


Итак, напомним еще раз назначение рассмотренных функций:

trigreduce
преобразует тригонометрическое выражение как сумма слагаемых, каждое из которых содержит один синус или косинус.
trigexpand
«раскрывает» аргументы тригонометрических функций, согласно правилам тригонометрических функций от суммы углов.
trigsimp
упрощает тригонометрическое выражение, используя тождество <m>sin^2x+cos^2x=1</m>.
trigsimp
мощная функция — упрощение рационального тригонометрического выражения.

Разложение на простые дроби: шаг за шагом[править | править код]

Ответы Maxima всегда показывают готовый результат (если, конечно, системе удалось его получить), но не путь его вычисления.

Например, разложение на простые дроби выполняется функцией «partfrac»:

partfrac ( 1/(x^2*(x^2 + 1)), x); → <m>{{1}\over{x^2}}-{{1}\over{x^2+1}}</m>

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

Итак, попробуем вычислить разложение в «полуавтоматическом» (или «полуручном») режиме.

В данном случае, разложение на простые дроби — это сумма следующих элементарных дробей:

p1: a/x;                  →   <m>{{a}\over{x}}</m>
p2: b/x^2;                   →   <m>{{b}\over{x^2}}</m>
p3: (c*x + d)/(x^2 + 1);        →   <m>{{c\,x+d}\over{x^2+1}}</m>

т. е.

p1 + p2 + p3;        →  <m>{{c\,x+d}\over{x^2+1}}+{{a}\over{x}}+{{b}\over{x^2}}</m>

Приводим к общему знаменателю, и поименуем числитель этой дроби, как «n»:

ratsimp(p1 + p2 + p3); →  <m>{{\left(c+a\right)\,x^3+\left(d+b\right)\,x^2+a\,x+b}\over{x^4+x^2
}}</m>
n: num(%);             →  <m>\left(c+a\right)\,x^3+\left(d+b\right)\,x^2+a\,x+b</m>

Теперь осталось решить уравнение, «приравняв» числитель исходной дроби и числитель «n»:

 solve ([coeff (n, x, 3) = 0,
         coeff (n, x, 2) = 0,
         coeff (n, x, 1) = 0,      →     [[a=0,b=1,c=0,d=-1]]
         coeff (n, x, 0) = 1],
         [a, b, c, d]);
 

Наконец, подставляем полученные коэффициенты:

at(p1 + p2 + p3, first(%));       →     <m>{{1}\over{x^2}}-{{1}\over{x^2+1}}</m>

Напомним, мы использовали следующие функции:

num
возвращает числитель дроби.
coeff
возвращает коэффициент при переменной в заданной степени.
solve
решает систему линейных уравнений.

Вычисление пределов[править | править код]

Maxima может вычислять пределы для широкого спектра как рациональных, так и трансцедентных функций.

Для рациональных функций, часто интересно посмотреть поведение функции на бесконечности:

limit ((x^2 - 1)/(x^2 + 1), x, inf); →  <m>\lim_{x\rightarrow \infty }{{{x^2-1}\over{x^2+1}}} \rightarrow 1 </m>

Можно раскрыть и предел типа «0/0»:

limit ((x^2 + x - 6)/(x^4 + x^3 - 19*x^2 + 11*x + 30), x, 2); →  <m>\lim_{x\rightarrow 2}{{{x^2+x-6}\over{x^4+x^3-19\,x^2+11\,x+30}}} 

\rightarrow -{{5}\over{21}}</m>

Maxima умеет вычислять пределы «справа» и «слева», и честно сообщает, когда предел в заданной точке не определен:

limit (tan(x), x, %pi/2, plus);   →  minf      <m>\lim_{x\rightarrow {{\pi}\over{2}}+0}{\tan x} \rightarrow -\infty</m>
limit (tan(x), x, %pi/2, minus);  →  inf       <m>\lim_{x\rightarrow {{\pi}\over{2}}-0}{\tan x} \rightarrow \infty</m>
limit (tan(x), x, %pi/2); → und

Вычисление рядов — конечных и бесконечных[править | править код]

Конечные суммы считаются сразу:

sum(i, i, 1, 100); → 5050  

Но бесконечные ряды, по умолчанию не суммируются, выражение остается в исходном символьном виде:

sum1:sum(1/x^2, x, 1, inf); →  <m>\sum_{x=1}^{\infty }{{{1}\over{x^2}}}</m>
sum2:sum(1/x^2, x, 1, inf); →  <m>\sum_{x=1}^{\infty }{{{1}\over{x^3}}}</m>

Чтобы выполнить суммирование, нужно указать опцию «simpsum=true»:

sum(1/x^2, x, 1, inf), simpsum=true; →  <m>{{\pi^2}\over{6}}</m>

Впрочем, над «непросуммированными суммами» можно выполнять операции, только операторы суммирования не будут сокращаться:

sum1+sum2; →  <m>\sum_{x=1}^{\infty }{{{1}\over{x^2}}}+\sum_{x=1}^{\infty }{{{1
}\over{x^3}}}</m>

Для сокращения сумм, нужно использовать функцию «sumcontract»:

sc:sumcontract(sum1+sum2); →  <m>\sum_{x=1}^{\infty }{{{1}\over{x^2}}+{{1}\over{x^3}}}</m>

Версия 5.9.3. содержит баг 1517061

Далее, можно просуммировать, как мы уже делали выше:

ev (sc, simpsum=true); →  <m>\zeta\left(3\right)+{{\pi^2}\over{6}}</m>

Внешние ссылки[править | править код]


По крайней мере часть этого текста взята с ресурса http://lib.custis.ru/ под лицензией GDFL.Список авторов доступен на этом ресурсе в статье под тем же названием.