Системы компьютерной алгебры

Mathematica допускает использование чисел четырех видов: целые, рациональные, вещественные и комплексные. Все типы чисел могут содержать любое количество цифр. Чтобы число рассматривалось как вещественное, оно должно содержать точку в его записи, даже если дробная часть равна нулю, например, 2. или 2.0.

Для перевода числа, заданного в системе счисления с произвольным основанием, в десятичную используется конструкция Основание^^Число, а для обратного перевода числа a из десятичной системы в систему с основанием n (где n не превышает 32) - функция BaseForm[a, n].

Инициализация переменных осуществляется при помощи операции =, для аннулирования значения переменной следует после знака равно указать символ . (точка).

При выполнении вычислений особая роль отводится символу % - он означает результат предыдущей операции. Комбинация символов %% соответствует результату операции, выполненной перед предыдущей, и так далее.

Для того чтобы "заставить" систему упростить выражение, используется функция Simplify. Ниже приведены примеры использования этой функции.

Функция Expand раскрывает скобки в выражении. Например, в результате выполнения команды Expand[(a + b)3] получится a3 + 3a2b + 3ab2 + b3.

К сожалению, функция Simplify не всегда выдает самый простой результат. В этом случае можно использовать функцию FullSimplify.

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

Обратите внимание на символ ; (точка с запятой) в конце ввода многочлена. Этот символ препятствует выводу на экран результата обработки программой Mathematica введенного выражения.

Разложение на множители, если это возможно, осуществляет функция Factor. Эта функция может работать и с тригонометрическими выражениями, но в этом случае нужно использовать дополнительную опцию Trig -> True:

In[10]:= Factor[x^4+8x^3+17x^2+16x+30]
Out[10]= (3 + x)(5 + x)(2 + x^2)

In[11]:= Factor[Sin[6x]/Sin[2x] + Cos[6x - Pi]/Cos[2x], 
               Trig -> True]
Out[11]= 2

Напомним, это многочленом P(x)степени n от переменной x называется выражение вида

a0xn+a1xn-1+...+ an-1x+an.

Для вынесения общего числового множителя в многочлене за скобки предназначена функция FactorTerms[poly, x], где poly есть многочлен от переменной x. Для получения списка коэффициентов при степенях x, начиная с нулевой, используется функция CoefficientList[poly, x].

In[12]:= FactorTerms[1232x^4+168x+144, x]
Out[12]= (8 ((18 + 21x + 154x^4))

In[13]:= CoefficientList[12x^4+68x+44,x]
Out[13]= {44,68,0,0,12}

Некоторые другие функции для работы с многочленами приведены в следующей таблице.

Пример

Пусть P1(x)= x4+2x3-4x2-5x-6. Определим, является ли число 2 корнем уравнения P1(x)=0. Известно, что многочлен делится без остатка на выражение x-x0, где x0 - корень уравнения. Найдем остаток от деления P1(x) на x-2:

PolynomialRemainder[x^4+2x^3-4x^2-5x-6, x-2 , x].

Результат равен 0, следовательно 2 - корень данного уравнения.

Задания

Вычислитe 2^-10 с точностью 20 знаков после запятой.

Упростите выражение

Разложите на множители выражение

x6-18x5+135x4-540x3+ 1215x2-1458x+729

.

Найдите остаток от деления многочлена P1(x) на x-1

Создание графических изображений

В системе Mathematica для построения графиков вида y=f(x) используется функция Plot, первым параметром которой является функция или список функций, графики которых должны быть построены, а вторым - список из трех величин: аргумент функций и концы промежутка, на котором должен быть построен график, например,

Plot[{Sin[x/2], 1.5 Sin[2x]}, {x, -3Pi,3Pi}]

Полученные изображения можно сохранить в форматах EPS, GIF, PDF, TIFF и некоторых других (смотри Help) при помощи функции Display. Ее первым аргументом является имя создаваемого файла, заключенное в кавычки, вторым - графическое изображение, а третьим - формат выходного файла, например,

Display["g1.gif", 
   Plot[ArcTan[x]-0.5, {x, -10, 10}],
   "GIF"]

Для каждого графика можно задать набор параметров PlotStyle, задающий цвет, толщину и тип линии графика. Опция RGBColor[r,g,b] (каждая составляющая есть число в диапазоне от 0 до 1) задает цвет графика в цветовой модели RGB. Для применения цветовой модели HSB (H - hue, тон; S - saturation, насыщенность; B - brightness, яркость) указывается опция Hue[h, s, b], (h, s и b из интервала [0; 1]). Возможно указание сокращенной формы этой опции - Hue[h], соответствующей значению Hue[h, 1, 1]. В следующем примере используются две различные цветовые модели для изображения двух графиков:

Plot[{Sin[x], Sin[2x]}, {x, -2Pi,2Pi},  
     PlotStyle -> {RGBColor[1,0,0], Hue[.4]}]

Другими опциями PlotStyle являются толщина линии (Thickness) и характеристики пунктирной линии (Dashing):

Plot[{x Sin[2x], x^3 - 3x +1}, {x, -4, 4},
     PlotStyle -> {{Hue[0], Thickness[.01]}, 
                  {Hue[.6], Dashing[{.02}]}}]

Функция Plot допускает задание множества других опций, относящихся как ко всему изображению в целом, так и к каждому из его составляющих в отдельности. Их полный список можно получить при помощи команды Options[Plot]. Вот примеры использования некоторых из них:

Пример

На следующем рисунке изображен график функции arctg x и две его горизонтальные ассимптоты (PI/2 и -PI/2). Опция PlotStyle задает параметры каждого графика в отдельности, а остальные опции задают внешний вид всего изображения в целом. Переменная p2 используется затем для сохранения изображения.

p2=Plot[{ArcTan[x], Pi/2, - Pi/2}, {x, -10, 10},
     PlotStyle -> {{RGBColor[0, 1, 0]}, 
                  {RGBColor[1, 0, 1], Dashing[{.015}]},
                  {RGBColor[1, 0, 1], Dashing[{.015}]}},
     Ticks -> {Automatic, {-Pi/2, -Pi/4, Pi/4, Pi/2}}, 
     PlotLabel->"f(x) = arctg(x)", 
     PlotRange ->{-2, 2}]

Для сохранения в файле g2.gif полученного изображения в формате GIF использовалась функция Display:

Display["g2.gif", p2, "GIF"]

Функция ParametricPlot используется для рисования графиков функций, заданных параметрически. Отметим, что Mathematica по умолчанию отображает графические объекты таким образом, что отношение высоты рисунка к его ширине равно величине, обратной золотому сечению (около 0.618034). Для сохранения естественных пропорций изображения следует указывать опцию AspectRatio -> Automatic. После этого круг выглядит действительно кругом, а не овалом. Ниже приведены примеры, демонстрирующие построение параметрически заданных кривых.

ParametricPlot[{Sin[t], Cos[t]}, {t, 0, 2Pi}]
ParametricPlot[{Sin[t], Cos[t]},
          {t, 0, 2Pi}, AspectRatio -> Automatic]
ParametricPlot[{Sin[3t], Cos[2t]},
          {t, 2Pi/3, 4Pi/3}, AspectRatio -> Automatic]
ParametricPlot[{{Sin[2t], Sin[3t]},
          {0.5*Sin[t], 0.5*Cos[t]}},
          {t, 0, 2Pi}, AspectRatio -> Automatic]

Для построения графиков в полярной системе координат следует предварительно загрузить модуль расширения Graphics, после чего использовать функцию PolarPlot:

<<Graphics`Graphics`            
PolarPlot[phi, {phi,0,10Pi}] 
PolarPlot[Cos[t], {t,0,2 Pi}]   
PolarPlot[Cos[a]^2 -Sin[a]^2,{a,0,2Pi}]
PolarPlot[Sin[5t],{t,0,2 Pi}]   
PolarPlot[1-Sin[fi],{fi,0,2Pi}]

Функции Plot3D и ParametricPlot3D являются основными инструментами программы Mathematica для изображения трехмерных графических объектов. Эти функции кроме трех обязательных параметров, задающих функцию и диапазоны изменения переменных, используют множество дополнительных опций. Отметим среди них только PlotPoints, указывающую на сколько частей следует разбить интервалы изменения переменных. Ниже приведены два примера использования этих функций.

Plot3D[Sin[xy], {x, 0, 4}, {y, 0, 4}, PlotPoints -> 40, Mesh ->False, FaceGrids ->All, AxesLabel -> {"Length", "Width", "Height"}]
ParametricPlot3D[ {u Cos[u](4 + Cos[u+v]), u Sin[u](4 + Cos[u+v]), u Sin[u + v]}, {u, 0, 3Pi},{v, 0, 3Pi}, PlotPoints ->{60, 12}]

Следующие изображения построены при помощи модуля расширения Graphics

Построение икосаэдра <<Graphics`Polyhedra` Show[Polyhedron[Icosahedron], Boxed -> False]
Лист Мебиуса - пример односторонней поверхности <<Graphics`Shapes` Show[Graphics3D[ MoebiusStrip[2, 0.5, 150], Boxed ->False]]
Сфера с отверстием <<Graphics`ContourPlot3D` ContourPlot3D[ Cos[Sqrt[x^2+y^2+z^2]], {x, -2, 2}, {y, -1.5, 2}, {z, -2, 2}, Boxed -> False, PlotPoints ->{3,15}]

Задания

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

Решение уравнений, систем и неравенств

Mathematica включает в себя средства поиска корней различных типов уравнений и систем. Чаще всего для этих целей используется функция Solve[eqn, var], где eqn - уравнение или система, задаваемая списком уравнений, var - переменная или список переменных, которые требуется определить. Результат вычислений - список корней (возможно пустой).

Палитра Basic Calculations в разделе Algebra содержит несколько шаблонов для ввода функции Solve (пункт Solving Equations ). Обратите внимание, что левая и правые части уравнения соединяются символами ==, например,

Solve[3 x + 9 == 0, x]
Solve[x^2 + 2b*x + c == 0, x]
Solve[x^2 + 1 == 0, x]
Solve[{x==1+2a*y, y==9+2x}, 
      {x,y}]

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

Программа позволяет также решать логарифмические и показательные уравнения. Натуральные логарифмы задаются функцией Log, а логарифмы по произвольному основанию n - в виде Log[n, expr], где expr - аргумент логарифмической функции, например, Log[2, 1024]. Для Эйлеровой константы (основания натурального логарифма) используется обозначение E.

Пример

Следующий фрагмент демонстрирует решение уравнения

In[12]:= Solve[Log[Sqrt[x]]==Sqrt[Log[x]], x]
Out[12]= {{x -> 1}, {x -> E^4}}

Для решения неравенств предварительно следует подключить дополнительный модуль расширения Algebra:

<<Algebra`InequalitySolve`

Имя функции, решающей неравенства, совпадает с названием подпакета - InequalitySolve. Для задания неравенств используют знаки > (строго больше), >= (больше или равно), < (строго меньше) и <= (меньше или равно), например,

In[13]:= InequalitySolve[x^2-1<0,x]
Out[13]= -1<x<1

In[14]:= InequalitySolve[x^2/(x^2-1)>=0,x]
Out[14]= x<-1 || x==0 || x>1

In[15]:= InequalitySolve[(2x-1)/(4(x+1))<Cos[5Pi/3],x]
Out[15]= x<-1

При решении неравенств, отличных от дробно-рациональных, система выдает предупреждение о возможном некорректном решении:

In[16]:= InequalitySolve[2^(x^2)>16,x]
InequalitySolve::"npi": "A nonpolynomial equation
       or inequality encountered. The solution
       set may be incorrect."
Out[16]= x<-2 || x>2

В тех случаях, когда невозможно решить заданное уравнение аналитически, можно приближенно вычислить значение корня следующим образом. Сначала при помощи функции Plot строятся графики левой и правой частей уравнения и по рисунку находится первое приближение x0. Затем для уточнения значения корня применяется функция FindRoot. Если L(x) и R(x) - соответственно левая и правая части уравнения, x - искомый корень, а x0 - его приближенное значение, определенное с помощью графика, то функция для поиска корня имеет вид:

FindRoot[L(x)==R(x), {x,x0}]

Пример

Найдем приближенное решение уравнения ex=x2. Построим графики функций ex и x2, выполнив команду Plot[{Exp[x], x^2}, {x, -1, 1}].

Мы получим изображение графиков функций ex и x2 на одном чертеже для значений аргумента x, изменяющегося на отрезке [-1,1].

Глядя на графики, видим, что в качестве начального приближения может быть взято значение x0 = -0.5.

Для уточнения решения нашего уравнения воспользуемся командой FindRoot[Exp[x]==x^2, {x, -0.5}] и получим приближенное значение корня x = -0.703467. Если требуется большая точность (но не более 16 знаков после запятой), то можно воспользоваться функцией N. Так, выполнив команду N[FindRoot[Exp[x]==x^2, {x, -0.5}], 16], мы получим приближенное значение корня с точностью до 16 знаков после запятой: x = -0.7034674295409824.

Пример

Как правило, тригонометрические уравнения имеют бесконечно много решений, поэтому можно указывать различные начальные приближения и получать приближенные значения разных корней уравнения. Так, выполнив команду FindRoot[Sin[x^2]==0, {x, 2}], найдем x = 1.77245, а изменив начальное приближение с 2.0 на 3.0, с помощью команды FindRoot[Sin[x^2]==0, {x, 3}] получим приближенное значение другого корня уравнения, а именно x=3.06998.

Задания

Решите системы уравнений:

1.

2.

Решите неравенства:

1.

2.

Найдите приближенно наименьший положительный корень уравнения 1/x2=5 cos x.

Найдите с точностью 12 знаков после запятой все корни уравнения (1 - x)/(x4 + 1) = sin x, принадлежащие отрезку [-1,4].

Пределы и ряды

Для вычисления пределов используется функция Limit. Необязательный параметр Direction применяется для вычисления односторонних пределов, которые используют, в частности, при исследовании функции на непрерывность.

Внимание: Direction -> 1 является указанием для вычисления предела слева, а Direction -> -1 - для предела справа.

Константа можно использовать палитру Basic Input.

Не все пределы могут быть вычислены таким образом. Дополнительные возможности предоставляет пакет расширений Calculus, подпакет Limit которого переопределяет встроенную функцию Limit. Для загрузки используется команда

<<Calculus`Limit`

Пример

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

Limit[Cos[n!]/n, n -> Infinity]

Пример

Используем односторонние пределы для исследования функции на непрерывность. Напомним, что функция называется непрерывной в точке x, если она определена в ней и ее значение в этой точке совпадает с левым и правым пределами функции при стремлении к x. Точка x является точкой устранимого разрыва, если пределы слева и справа существуют (принимают конечные значения) и равны, но значение функции f(x) не существует или отлично от них. Точка x называется точкой разрыва первого рода, если пределы не совпадают, но конечны. Наконец, x есть точка разрыва второго рода, если в ней не существует хотя бы одного конечного предела.

Для функции sin(x)/x точка x = 0 является точкой устранимого разрыва, так как существуют пределы слева и справа, равные 1.

Довольно часто в процессе решения различных задач математического анализа приходится сталкиваться с рядами и разложением аналитических функций в ряды. Поэтому в программе Mathematica имеется обширный набор функций для решения подобного рода проблем. Для получения частичной суммы ряда используется функция Sum[f, {i, imin, imax}]. Здесь f - общий член ряда, а imin и imax - границы суммирования.

Большинство достаточно гладких функций довольно точно локально аппроксимируются рядом Тейлора. Для нахождения n -го многочлена Тейлора в окрестности точки x0 применяется функция Series[f, {x, x0, n}], а удаление остаточного члена ряда выполняется функциями Collect и Normal

Пример

Для того чтобы наглядно представить себе точность апроксимации функции многочленами Тейлора, построим на одном рисунке график функции sin(x)/x и соответствующие ей многочлены Тейлора степеней 3 и 6.

f=Sin[x]/x;
f1= Normal[Series[f, {x, 0, 3}]];
f2= Normal[Series[f, {x, 0, 6}]];
Plot[{f, f1, f2}, {x, -10, 10}, 
  PlotStyle -> {
     {Hue[0], Thickness[.01]},
     {Hue[.6], Dashing[{.03}]},
     {Hue[.9], Dashing[{.01}]}}]

Задания

Вычислите пределы:

Найдите односторонний предел

Исследуйте функции на непрерывность:

Дифференцирование и интегрирование

Для вычисления производной функции f по переменной x используется команда D[f, x]. Нахождениe n -ой производной выглядит так: D[f,{x,n}]. На панели BasicInput имеются кнопки с шаблонами для дифференцирования.

Пример

Построим графики функции arctg x - 0.5, ее первой и второй производных. Обратите внимание на определение функций в этом примере. Символ подчеркивания после имени переменной означает, что определяется функция от указанной переменной. Процедура Clear уничтожает данное ранее определение.

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

In[1]:= Integrate[x^2/(4x^6+1), x]
Out[1]:= (1/6) ArcTan[2x^3]

In[2]:= Integrate[1/Sqrt[(4-x)^3],{x,-Infinity,0}]
Out[2]= 1

Зачастую полезно применить функцию Simplify к результатам вычисления, например,

Задание

Операции с матрицами

Матрица задается как список списков элементов строк. Для представления списка в традиционной матричной форме используется функция MatrixForm. Процедура Det позволяет вычислять определитель. С помощью функции Dot находится произведение матриц и векторов, но можно использовать и символ . (точка). Inverse находит обратную матрицу, а Transpose используется для транспонирования матрицы (превращения каждой строки в столбец с тем же номером).

Для решения систем линейных уравнений можно использовать функцию LinearSolve или матричный метод, при котором матрица X, содержащая значения неизвестных, находится по формуле X = A-1x B, где x означает операцию умножения матриц.

Пример

Используем указанные способы при решении следующей системы уравнений:

Итак, решением данной системы является тройка {-1, 2.5, 1.5}.

Задания

Найдите произведение матриц ,

Транспонируйте матрицу B и найдите ее определитель.

Логические выражения

Программа !, логическое отрицание). Импликация задается функцией из палитры BasicInput. Функция LogicalExpand осуществляет упрощение логических выражений.

Пример

Определим значение логического выражения !(x>z) !(x=y), если x=3, y=5, z=2

Задание

Получите законы дополнения ( !(!a) = a ), идемпотентности ( a a = a; a || a = a ) и поглощения ( a (a || b) = a ) при помощи функции LogicalExpand.

Система символьной математики Maxima

Maxima - еще одна программа для выполнения математических вычислений, символьных преобразований, а также построения разнообразных графиков. Сложные вычисления оформляются в виде отдельных процедур, которые затем могут быть использованы при решении других задач. Система Maxima распространяется под лицензией GPL и доступна как пользователям ОС Linux, так и пользователям MS Windows.

Для работы с данной системой в ОС Linux следует в окне shell набрать команду maxima или xmaxima для запуска ее графической оболочки. Другим удобным инструментом для работы с системой Maxima является программа texmacs. На панели инструментов этой программы располагается кнопка с изображением монитора, нажатиe на которую открывает меню выбора интерактивной сессии. Выбор пункта maxima позволит начать сеанс работы с этой программой.

При отображении результатов вычислений эта оболочка использует стандартные математические обозначения, в то время как xmaxima или maxima - только символы из таблицы ASCII-кодов.

При старте выводится некоторая информация о системе и "метка" ( C1 ). Каждый ввод и вывод помечаются системой и затем могут быть использованы снова. Символ C (от command) используется для обозначения команд, введенных пользователем, а D (от display) - при выводе результатов вычислений.

Для инициализации процесса вычислений следует ввести команду, затем символ ; (точка с запятой) и нажать клавишу Enter. Если не требуется вывод полученной информации на экран, то вместо точки с запятой используется символ $. Обратиться к результату последней команды можно с помощью символа %. Для повтора ранее введенной команды, скажем ( C2), достаточно ввести два апострофа и затем метку требуемой команды, например, "C2.

Система Maxima не обращает внимание на регистр введенных символов в именах встроенных констант и фунций. Запись sin(x) эквивалентна записи SIN(x), но при выводе результатов в текстовом режиме используются заглавные буквы. Регистр букв, однако, важен при использовании переменных, например, Maxima считает x и X разными переменными.

Для стандартных математических констант используются следующие обозначения: %e (или %E ) для основания натуральных логарифмов, %i ( %I ) для мнимой единицы (квадратный корень из числа -1) и %pi ( %PI ) для числа $$\pi$$.

Присваивание значения какой-либо переменной осуществляется с помощью знака : (двоеточие), а символ = (равно) используется при задании уравнений или подстановок.

Функция kill аннулирует присвоенные ранее значения переменных. Параметр all этой функции приводит к удалению значения всех переменных, включая метки Ci и Di.

Для завершения работы с системой применяется функция quit();, а прерывание процесса вычислений осуществляется путем нажатия комбинации клавиш Ctrl+c (после чего следует ввести :q для возврата в обычный режим работы).

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

Работа с выражениями

При записи математических выражений могут использоваться четыре стандартные арифметрические операции (+, -, *, /) и операция возведения в степень (^, ^^ или **). Приоритет этих операций традиционен, для изменения порядка вычислений следует использовать круглые скобки. Кроме чисел выражения могут содержать результаты вычислений математических функций. Аргументы функций указываются в круглых скобках, например, запись sqrt(5) означает корень квадратный из числа 5. Если в результате расчета получается дробное выражение, то оно и выводится в виде обыкновенной дроби. Иррациональные числа, входящие в выражение, представляются в символьном виде.

Для раскрытия скобок используется функция expand. Команда ev позволит получить численное значение выражения. Ее первый аргумент есть вычисляемое выражение, а второй - опция numer. Напомним, что символ % означает результат предыдущего вычисления.

На рисунке изображен пример работы с оболочкой texmacs. Дальнейший вывод соответствует работе с программой maxima ( xmaxima ).

Допускается более удобная форма функции ev, требующая указания только ее аргументов:

C5) 29*sqrt(2) + 41, numer;
(D5)                      82.01219330881976

По умолчанию результат содержит 16 значащих цифр. Для вывода числа в экспоненциальной форме используется функция bfloat:

(C6) bfloat(d3);
(D6)                       8.201219330881976B1

Запись mBn есть сокращенная форма выражения m*10ⁿ.

Количество значащих цифр в представлении числа определяется специальной переменной FPPREC. Увеличение ее значения приводит к возрастанию точности результата, например,

(C7) fpprec;
(D7)                                  16

(C8) fpprec:100;
(D8)                                  100

(C9) "'c5;
(D9) 8.20121933088197607657604923686248525#
030775305167218616484631047782707024434954#
8350683851114422615155B1

Символ # в конце выводимой строки означает, что число не уместилось на одной строке и его оставшаяся часть переносится на следующую. В последнем примере мы использовали повторение ранее введенной команды ( "c5 ).

Система Maxima может работать с числами произвольной длины и точности:

(C10) 100!;
(D10)  933262154439441526816992388562667004#
9071596826438162146859296389521759999322991#
5608941463976156518286253697920827223758251#
185210916864000000000000000000000000

Пример

Вычислим число PI с точностью 200 знаков после запятой:

(C11) %PI, numer;
(D11) 		 3.141592653589793

(C12) fpprec:200;
(D12) 		        200

(C13) bfloat(%PI);
(D13) 3.141592653589793238462643383279502884#
19716939937510582097494459230781640628620899#
86280348253421170679821480865132823066470938#
44609550582231725359408128481117450284102701#
9385211055596446229489549303819B0

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

(C14) FPPREC:16;
(D14) 	                              16

(C15) bfloat(%PI);
(D15) 	             3.141592653589793B0

(C16) block([FPPREC:100], bfloat(%PI));
(D16) 3.14159265358979323846264338327950#
2884197169399375105820974944592307816406#
286208998628034825342117068B0

(C17) "c15;
(D17)                3.141592653589793B0

В Maxima доступны прямые и обратные тригонометрические функции: sin (синус), cos (косинус), tan (тангенс), cot (котангенс), asin (арксинус), acos (арккосинус), atan (арктангенс), acot (арккотангенс). Кроме них имеются две менее известные функции - секонс ( sec x = 1/cos x ) и косеконс ( csc x = 1/sin x ).

(C18) block([FPPREC:100],sin(bfloat(%PI)));
(D18) - 2.570579198029723711689569064713781#
2738478411385601247337570449378000209830368#
31739403591933813645377B-101

(C19) block([FPPREC:100],sin(%PI)) ;
(D19)                                 0

(C20) sin(%PI/6)^2 +cos(%PI/6)^2;
(D20)                                 1

(C21) sec(%PI/3);
(D21) 		                      2

(C22) asin(1/2);
                                   %PI
(D22)                              ---
                                    6

К сожалению, Maxima не может в символьном виде вычислить значения обратных тригонометрических функций в некоторых точках:

C23) asin(sqrt(3)/2);
			       SQRT(3)
(D23) 			  ASIN(-------)
				  2
(C24) %, expand;
			       SQRT(3)
(D24)			  ASIN(-------)
				  2

Для вычисления натурального логарифма используется функция log:

(C25) log(%E^2);
(D25) 				 2
(C26) 5*log(a)+6*log(b);
(D26) 			 6 LOG(b) + 5 LOG(a)
(C27) logcontract(%);
				   5  6
(D27) 			      LOG(a  b )

Алгебраические преобразования

Особо важную роль играет способность системы Maxima производить разнообразные символьные преобразования алгебраических выражений. Следующий пример демонстрирует раскрытие скобок в них:

C1) p1:x^2-1;
				 2
(D1) 			        x  - 1

(C2) p2:x-1;
(D2) 				 x - 1

(C3) expand(p1*p2);
			     3	  2
(D3) 			    x  - x  - x + 1

(C4) expand((p1+p2)^2);
		       4      3	     2
(D4) 		      x  + 2 x  - 3 x  - 4 x + 4

Функцию divide можно использовать для нахождения частного и остатка от деления одного многочлена на другой:

(C5) divide(p1*p2, p1);
(D5) 			      [x - 1, 0]

Функция gcd определяет наибольший общий делитель многочленов, а factor осуществляет разложение многочлена на множители:

(C6) gcd(x^3-1, x^2-1, x-1);
(D6)                         x - 1

(C7) factor(x^8-1);
        		     2	      4
(D7)	   (x - 1) (x + 1) (x  + 1) (x  + 1)

Подстановка какого-либо выражения вместо переменной осуществляется при помощи операции =. Например, заменим все вхождения x в выражении на 5/z:

(C8)  x^4+3*x^3-2*x, x=5/z;
			 10   375   625
(D8) 		       - -- + --- + ---
			 z     3     4
			      z	    z

Функция ratsimp выносит за скобки наибольший общий делитель:

C9) ratsimp(%);
			   3
		       10 z  - 375 z - 625
(D9) 		     - -------------------
				4
			       z

Используя функцию assume (to assume - допускать), можно при вычислениях учитывать дополнительные условия, задаваемые неравенствами:

(C10) sqrt(x^2);
(D10) 			       ABS(x)
(C11) assume(x<0);
(D11) 			       [x < 0]
(C12) sqrt(x^2);
(D12) 				 - x

Функция forget (to forget - забывать) снимает все ограничения, наложенные при помощи assume:

C13) forget(x<0);
(D13) 			       [x< 0]
(C14) sqrt(x^2);
(D14) 			       ABS(x)

Maxima легко оперирует тригонометрическими выражениями. Так, функция trigexpand использует формулы преобразования сумм двух углов для представления введенного выражения в как можно более простом виде:

(C15) sin(u+v)*cos(u)^3;

		   3
(D15) 		COS (u) SIN(v + u)

(C16) trigexpand(%);

	   3
(D16) 	COS (u) (COS(u) SIN(v) + SIN(u) COS(v))

Функция trigreduce преобразует тригонометрическое выражение к сумме элементов, каждый из которых содержит только единственный sin или cos:

C17) trigreduce(%);
 
            SIN(v + 4 u) + SIN(v - 2 u)   
(D17)       --------------------------- 
                         8                
 
            3 SIN(v + 2 u) + 3 SIN(v)
          + -------------------------
                         8

Функции realpart и imagpart возвращают действительную и мнимую часть комплексного выражения:

(C18) z1:-3+%i*4;
(D18) 		        4 %I - 3

(C19) z2:4-2*%i;
(D19) 		        4 - 2 %I

(C20) z1*z2;
(D20) 	   (4 - 2 %I) (4 %I - 3)

(C21) expand(%);
(D21) 		       22 %I - 4

(C22) realpart(''c20);
(D22) 		             - 4

(C23) imagpart("c20);
(D23) 		              22

Решение уравнений

Maxima может решать уравнения и системы алгебраических уравнений с помощью функции solve. Равная нулю правая часть уравнения может быть опущена:

(C1) solve (x^2=1, x);
(D1) 			        [x = - 1, x = 1]

(C2) solve(x^2-1,x);
(D2) 			        [x = - 1, x = 1]

(C3) solve(log(x+3)=1, x);
(D3) 				  [x = %E - 3]

При решении тригонометрических уравнений выдается только одно из бесконечного множества возможных решений:

(C4) solve(sin(x)-1, x);
SOLVE is using arc-trig functions to get
a solution. Some solutions will be lost.
				 %PI
(D4) 			    [x = ---]
				  2

В следующем примере функция solve используется для решения системы из трех уравнений с тремя неизвестными:

(C5) s:[x+y+z=3, x+2*y-z=2, x+y*z+z*x=3];
(D5) 	      [z + y + x = 3, - z + 2 y + x = 2,
               y z + x z + x = 3]
(C6) solve(s, [x,y,z]);
(D6) 	       [[x = 1, y = 1, z = 1], 
                [x = 7, y = - 3, z = - 1]]

Если уравнение не имеет решений на множестве действительных чисел, то Maxima ищет решения среди комплексных чисел:

(C7) solve(x^2+1,x);
(D7) 			       [x = - %I, x = %I]

Задание

Решите уравнение ln(tg x) = 0.

Определение функций

Для определения функции следует указать ее имя, аргумент (или аргументы), заключенный в круглые скобки, символы := и саму функцию, например,

C1) p(x):= x^2+3*x-8;
			             2
(D1) 		            p(x) := x  + 3 x - 8
(C2) p(2*x);
			             2
(D2) 		                  4 x  + 6 x - 8
(C3) factor(p(2*x-4));
			            2
(D3) 		              2 (2 x  - 5 x - 2)
(C4) f(x):=sin(x);
(D4) 		                  f(x) := SIN(x)
(C5) g(x):=log(x);
(D5) 		                  g(x) := LOG(x)
(C6) f(g(x));
(D6)		                     SIN(LOG(x))
(C7) g(f(x));
(D7)		                     LOG(SIN(x))
(C8) dist(x,y):=sqrt(x^2+y^2);
	                                 2    2
(D8)	             dist(x, y) := SQRT(x  + y )
(C9) dist(3,4);
(D9) 	                                       5

Функция может быть задана и подобным образом:

(C10) eta(x,y):= if x = y then x 
               else (if x > y then x-y else x+y);
(D10) 	    eta(x, y) := IF x = y THEN x ELSE 
              (IF x > y THEN x - y ELSE x + y)
(C11) eta(5,6);
(D11) 				  11
(C12) eta(eta(7,7), eta(1,2));
(D12) 				  4

Построение графиков

Графики строятся при помощи функции plot2d, например,

(C1) plot2d([2*sin(x), cos(x)], [x, -2*%PI, 2*%PI]);

Для построения изображений трехмерных объектов используется функция plot3d. Вот как, например, запишется функция для рисования листа Мебиуса:

(C2) plot3d([cos(x)*(3+y*cos(x/2)), 
  sin(x)*(3+y*cos(x/2)), y*sin(x/2)],
  [x,-%pi,%pi],[y,-1,1],['grid,40,15]);

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

(C3) Plot3d(r^.33*cos(th/3),[r,0,1],[th,0,6*%pi],
     ['grid,12,80],
     ['transform_xy,polar_to_xy],
     ['view_direction,1,1,1.4],
     ['colour_z,true]);

Операции математического анализа

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

Для нахождения производной используется функция diff, первым аргументом которой является функция, вторым - переменная, по которой производится дифференцирование, и третьим (необязательным) - порядок производной:

(C1) f:(x-2*sqrt(x))/x^2;
                         x - 2 SQRT(x)
(D1) 	                 -------------
		               2
		              x

(C2) diff(f, x);
           1
      1 - -------
	  SQRT(x)   2 (x - 2 SQRT(x))
(D2)  ----------- - -----------------
	   2		    3
	  x		   x

(C3) expand("c2);
		             3     1
(D3) 		            ---- - --
		             5/2    2
		            x      x

(C4) g:x^6;
		                   6
(D4) 		                  x

(C5) diff(g, x, 1);
		                   5
(D5) 		                6 x

(C6) diff(g, x, 4);
		                   2
(D6) 		              360 x

При вычислении кратных производных по нескольким переменным после указания функции перечисляются переменные дифференцирования с указанием соответствующих кратностей, например,

(C7) diff(x^6*y^3, x, 4, y, 2);
		                   2
(D7) 		          2160 x  y

Функция integrate позволяет вычислять интегралы. Для нахождения неопределенного интеграла после функции указывается единственный аргумент - переменная интегрирования:

(C8) f:x^2/(4*x^6+1);
		                2
		               x
(D8) 		            --------
	                       6
                            4 x  + 1

(C9) integrate(f, x);
			          3
         		  ATAN(2 x )
(D9) 	                  ----------
		              6

Maxima в случае неоднозначного ответа может задавать дополнительные вопросы, как в следующем примере:

(C10) integrate(x^n,x);
Is  n + 1  zero or nonzero?
nonzero;
                               n + 1
		              x
(D10) 		              ------
		              n + 1
(C11) integrate(x^n,x);
Is  n + 1  zero or nonzero?
zero;
(D11)                         LOG(x)

Можно использовать функцию assume для задания дополнительных условий (не забывайте затем удалить наложенные ограничения):

(C12) assume(notequal(n,-1));
(D12) 			 [NOT EQUAL(n, - 1)]
(C13) integrate(x^n,x);
			      	 n + 1
				x
(D13) 			        ------
				n + 1
(C14) forget(notequal(n,-1));
(D14) 			 [NOT EQUAL(n, - 1)]
(C15) integrate(x^n,x);
Is  n + 1  zero or nonzero?

zero;
(D15) 			        LOG(x)

Для нахождения определенного интеграла следует указать дополнительные агрументы - пределы интегрирования:

(C16) integrate(x^2, x, 0, 6);
(D16) 	                              72

(C17) integrate(sin(x), x, 0, %PI);
(D17) 				       2

(C18) integrate(integrate(x*y, x, 1, 3), y, 0, 4);
(D18) 				      32

Maxima допускает задание и бесконечных пределов интегрирования. Для обозначения бесконечности используется переменная INF ( inf ):

(C19) integrate(1/x^2, x, 1, inf);
(D19)                          1  

(C20) integrate(1/(1+x^2), x, -inf, inf);
(D20) 				      %PI
 
(C21) integrate(1/x, x, 0, inf);
 
Integral is divergent  - an error.  
Quitting. To debug this try DEBUGMODE(TRUE);)

В последнем примере система сообщила о невозможности вычисления интеграла, т. к. он расходится (is divergent).

При вычислении достаточно сложных интегралов ответ не всегда будет представлен в наиболее простом виде. В следующем примере Maxima не может в символьном виде получить ответ, равный PI/4:

C22) g:1/sqrt(2-x^2);
                		    1
(D22) 	                      ------------
			                2
                              SQRT(2 - x )

(C23) integrate(g,x, 0,1);

                                  SQRT(2)
(D23) 	                     ASIN(-------)
		                     2

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

(C24) sum(1/n^2,n,1,inf);
	                        INF
		               ====
		               \      1
(D24) 		                >    --
		               /      2
		               ====  n
	                       n = 1
(C25) %,simpsum;
			              2
		                   %PI
(D25) 		                   ----
		                    6

Maxima способна находить разложение функций в ряд Тейлора. Получим многочлен Тейлора порядка 4 для функции f(x)=ln x в точке x=1:

(C26) g:log(x);
(D26)                            LOG(x)
(C27) taylor(g,x,1,4);
		       2	  3	     4
		 (x - 1)    (x - 1)    (x - 1)
(D27)/T/ x - 1 - -------- + -------- - -------- + ...
		    2	       3	  4

Для вычисления пределов используется функция limit:

(C28) limit(1/x,x,inf);
(D28) 		                   0

Для вычисления односторонних пределов используется дополнительный параметр, принимающий значение plus для вычисления предела справа и minus - слева.

Пример

Исследуем на непрерывность функцию arctg(1/(x-4)). Эта функция не определена в точке x = 4. Вычислим пределы справа и слева:

(C28) limit(atan(1/(x-4)), x, 4, plus);
				      %PI
(D28) 				      ---
				       2
(C29) limit(atan(1/(x-4)), x, 4, minus);
				       %PI
(D29) 				     - ---
					2

Как видим, точка x = 4 является точкой разрыва I рода для данной функции, так как существуют пределы слева и справа, равные - PI/2 и PI/2 соответственно.

Задания

Матричные вычисления

Maxima позволяет легко манипулировать матрицами. В следующем примере задаются две матрицы, которые затем складываются (+) и перемножаются (.):

(C1) A:matrix([1,2],[3,4]);
				    [ 1	 2 ]
(D1) 				    [ 	   ]
				    [ 3	 4 ]
(C2) B:matrix([1,1],[1,1]);
				    [ 1	 1 ]
(D2) 				    [ 	   ]
				    [ 1	 1 ]
(C3) A + B;
				    [ 2	 3 ]
(D3) 				    [ 	   ]
				    [ 4	 5 ]
(C4) A . B;
				    [ 3	 3 ]
(D4) 				    [ 	   ]
				    [ 7	 7 ]

Функция determinant вычисляет определитель матрицы.

(C5) determinant(A);
(D5) 				  - 2

(C6) determinant(matrix([a,b],[c,d]));
(D6) 				  a d - b c

Транспонирование матрицы осуществляется функцией transpose.

(C7) transpose(A);
			       [ 1  3 ]
(D7) 			       [      ]
			       [ 2  4 ]

Для получению обратной матрицы используется операция ^^-1 или функция invert.

(C8) A^^-1;
			[ - 2    1  ]
	                [           ]
(D8) 			[  3      1 ]
			[  -    - - ]
			[  2	  2 ]
(C9) invert(A);
                        [ - 2    1  ]
                        [           ]
(D9)                    [  3      1 ]
                        [  -    - - ]
                        [  2      2 ]

Как известно, каждый элемент bij обратной матрицы B = A-1 получается делением алгебраического дополнения Aij соответствующего элемента исходной матрицы на ее определитель |A|. Для того чтобы вынести 1/|A| в качестве сомножителя применяется функция detout.

(C10) invert(A), detout;
			      [	 4   - 2 ]
			      [ 	 ]
			      [ - 3   1	 ]
(D10) 			    - ------------
				   2

Убедимся в правильности полученного результата, умножив A на обратную к ней матрицу:

(C11) A . d9;
			       [ 1  0 ]
(D11) 			       [      ]
			       [ 0  1 ]

Будьте внимательны: в результате выполнения операции ^-1 получится матрица, каждый элемент которой обратен элементу исходной, а не обратная матрица.

(C12) A^-1;
				    [	 1 ]
				    [ 1	 - ]
				    [	 2 ]
(D12) 				    [ 	   ]
				    [ 1	 1 ]
				    [ -	 - ]
				    [ 3	 4 ]

Использование матриц позволяет легко решать системы линейных уравнений с несколькими переменными. Пусть . означает матричное умножение.

Пример

Решим следующую систему уравнений матричным способом.

Сначала заполним соответствующие матрицы, а затем получим матрицу результатов:

(C14) A:matrix([1, 2, 1], [2, 1, 1], [1, 3, 1]);

			  [ 1  2  1 ]
			  [	    ]
(D14) 			  [ 2  1  1 ]
			  [	    ]
			  [ 1  3  1 ]

(C15) B:matrix([0, 1, 0]);
(D15) 			  [ 0  1  0 ]

(C16) (A^^-1).B;
			    [  1  ]
			    [	  ]
(D16) 			    [  0  ]
			    [	  ]
			    [ - 1 ]

Задания

Вычислите произведение матриц

Найдите определители матриц

Для матрицы D найдите обратную, после чего проверьте, что в результате их произведения получается единичная матрица.

Решите следующую систему уравнений матричным способом

Программирование в системе Maxima

До сих пор мы использовали систему Maxima в интерактивном режиме, подобно калькулятору. Если часто приходится выполнять определенную последовательность вычислений, то лучше оформить ее в виде программы, которая затем вызывается в случае надобности. Ниже приводится небольшая программа для нахождения критических точек функции f(x). Пользователю предлагается ввести функцию f, после чего вычисляется производная введенной функции и при помощи функции solve решается уравнение fx = 0. Программа записывается в текстовый файл и затем загружается в систему Maxima при помощи функции batch. Приведем текст программы:

/* -----------------------------------------------------
   Это файл critpoints.max: 
   (Maxima использует комментарии в стиле C)
   -------------------------------------------------- */

critpts():=(
   print("Программа нахождения критических точек"),

/* Запрос на ввод функции     */
   print("Введите функцию f(x):"),
   f:read(),

/* Печать введенной функции (для контроля) */
   print("f = ",f),

/* В переменную eq помещаем значение производной */
   eq:diff(f,x),

/* Решаем уравнение */
   solve(eq, x)   
)$

Программа состоит из единственной функции (без аргументов), которая называется critpts. Команды отделяются друг от друга запятыми. Вот пример выполнения программы:

(C1) batch("critpoints.max");

batching #p/home/test/critpoints.max
(C2) critpts() := (PRINT("Программа #
нахождения критических точек"), 
PRINT("Введите функцию f(x):"), 
f : READ(), 
PRINT("f = ", f), 
eq : DIFF(f, x), 
SOLVE(eq, x))

(C3) critpts() ;
Программа нахождения критических точек
Введите функцию f(x): 
(x+2)/(x^2+1);

                  x + 2
            f =  ------ 
                  2
                 x  + 1

(D3) [x = - SQRT(5) - 2, x = SQRT(5) - 2]