Дифференцирование в Maple "своими силами"

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

 

1. (a+b)'=a'+b' - дифференцирование суммы,
2. (a b)'=a' b+a b' - дифференцирование произведения,
3. (a / b)'=(a' b-a b') / b^2 - дифференцирование отношения,
4. f(g(x))'=f ' (g(x)) g'(x) - дифференцирование сложной функции,

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

Фактически, настоящей проблемой при дифференцировании является упрощение результата, так как если его не упрощать, то производная от 2x+1 формально выдается в виде 0x+2*1+0. В силу этого, в современных универсальных программах алгоритм дифференцирования стал значительно сложнее, так как в данных пакетах используется достаточно широкий класс кусочно-непрерывных функций.
Алгоритм дифференцирования достаточно прост и, поэтому, попробуем реализовать этот алгоритм для дифференцирования полиномов от одной переменной

 

P = a_nxⁿ + a_n-1x^(n-1) + ... + a₀

Замечание: Работа компьютера напоминает работу мельницы: засыпаешь зерно - получаешь муку, а засыпаешь камни - получаешь песок. Это особенно важно при работе с системами символьных вычислений. Перед началом работы с Maple необходимо понять цель работы, продумать различные подходы к решению проблемы и,как это не пародоксально, скорее всего уже знать или предполагать конечный ответ.

 

Итак, для того чтобы продифференцировать полином нам надо:

1. отличать в выражении P переменную x от коэффициентов a_k и показателей степеней n;
2. различать математические операции - сумму, произведение и степень;
3. задать правила дифференцирования. 

 

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

Пример 1 ( команды op - nop).

v := f(x,y,z);

nops(v);

op(0,v);

op(2,v);

Более сложный пример

w := f(g(a,b),h(c,d));

op(1,op(2,w));

op([2,1],w);

op([-1,-1],w);

Таким образом можно выделять более сложные части выражений.

s := series(sin(x),x=2,6);

op(0,s);

В последнем примере мы раздожили функцию sin(x) в окрестности точки x=2 с точностью O(x⁶).

 

Далее нам надо научиться определять тип переменных и выражений. Для этого используется логическая команда type(e,t) где e - выражение, а t - имя типа выражения. В ответе получаем либо "правда", либо "ложь". Типов выражений очень много и из простых типов можно составлять составные типы. Заметим, что команда whattype(e), которая должна выдавать тип выражения значительно уступает команде type по своим возможностям. Рассмотрим несколько примеров.

 

Пример 2 ( типы выражений и переменных).

type( 2, integer );

type( a + b, `+` );

type( a * b, `+` );

type( a and b, `and` );

Это простые типы, рассмотрим составные

type(x^(-2),name^integer);

type(x^(-2),name^posint);

type(x^(-2),algebraic^integer);

type(x+y^(1/2)+1,`&+`(name,radical,integer));

type(a*b*c,`&*`(algebraic,algebraic,algebraic));

type(exp(x),exp(name));

type([x,1],[name,integer]);

 

Так как type булевская команда, то изучим сразу же условный оператор if вида

if условие a then выполняемое выражение A

elif условие b then выполняемое выражение B

else выполняемое выражение C

fi;

 

Эта конструкция дает возможность в зависимости от условий a и b выполнять выражения A, B или C, в качестве которых может выступать любая последовательность Maple-команд.

 

Пример 3 (условный оператор).

a := 3; b := 5;

if (a > b) then a else b fi;

Более сложная конструкция

5*(Pi + `if`(a > b,a,b));

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

 

Пример 4 (процедура дифференцирования полиномов).

В заголовке процедуры мы присваеваем ей имя df. Входными параметрами будут имя дифференцируемого выражения p и имя независимой переменной x по которой мы его дифференцируем. После заголовка следует описание прцедуры, отделенное пробелом.

 

df:=proc(p:algebraic,x:name)

      local u,v;

      if type(p,numeric) then 0

      elif type(p,name) then

      if p=x then 1 else 0 fi

      elif type(p,`+`) then map (df,p,x)

      elif type(p,`*`) then

           u:=op(1,p); v:=p/u;

           df(u,x)*v+df(v,x)*u;

      elif type(p,anything^integer) then

           u:=op(1,p); v:=op(2,p);

           v*df(u,x)*u^(v-1);

     else ERROR(`Ошибка при дифференцировании полинома`,p);

     fi;

end:

 

Разберем структуру процедуры.

В первой строке мы определяем имя процедуры, типы входных данных p и x и вводим две внутренних переменных u,v (это можно сделать и по умолчанию).

Во второй строке задано правило: производная от числа (тип numeric) есть ноль.

В третьей строке задано правило: производная от переменной (тип name) есть dx / dx =1 или dy / dx=0.

В четвертой строке задано правило 1 о дифференцировании суммы (тип `+`). Для дифференцирования слагаемых применяется таже процедура df, к которой мы обращаемся реккурсивно, используя команду map (данная команда позволяет применить процедуру df ко всем слагаемым сразу).

В пятой строке задано правило 2 о дифференцировании произведения (тип `*`). Для дифференцирования множителей, которым присваиваются имена u и v, применяется процедура df . (в отличие от предыдущей строки мы вынуждены делать это явно).

В шестой строке задано правило дифференцирования степени переменной (тип anything^integer), при этом мы так же явно выделяем части выражения.

В седьмой строке для аварийного выхода из процедуры в случае возникновения ошибки используется команда ERROR (именно из больших букв).

 

Теперь можно продифференцировать что нибудь

g:=y^5+x*y^4+z*y^3;

df(g,y);

 

Встроенные процедуры дифференцирования

 

Для дифференцирования выражения e по переменным x, y и т.д. используют следующие две команды в среде Maple:

Diff, diff, - как обычно, команды с большой и малой буквы отвечают символу и значению т.е. процедуре отложенного и немедленного выполнения;

D - позволяет найти дифференциал функции и дифференцировать процедуры;

implicitdiff - позволяет дифференцировать функции, заданные уравнением.

 

Пример 5 (дифференцирование )

Diff(sin(x),x);

diff(sin(x),x);

Diff(sin(x),x$3);

value(");

В последнем примере мы задали порядок дифференцирования.

Если функция не задана явно, то

diff(f(x),x);

diff(f(x,y),x,y);

Аналогично вычисляются и дифференциалы

D(sin);

D(ln);

D(f);

D(f)(0); # Производная от f вычисляется в точке 0

D(D(f));

(D@@2)(f); # Дифференциал старшего порядка как композиция простых дифференциалов

D(f@g); # Дифференциал от композиции двух функций

Эта комманда позволяет дифференцировать и составные выражения. Например, предположим, что a -это число

assume(a,constant);

D(2*sin(a)*f);

Дифференциал можно преобразовать в производную и обратно

f := D(y)(x) - a*D(z)(x);

convert(f, diff, x);

f := diff(y(x), x$2);

convert(f, D);

Рассмотрим уравнение, которое определяет x(y) или y(x)

f := x^2+y^3=1

Продифференцируем y по x или x по y

implicitdiff(f,y,x);

implicitdiff(f,x,y);

implicitdiff(f,y,z);

implicitdiff(f,y(x),x);

implicitdiff(f,y,x,x);

 

При желании можно посмотреть как реализованы процедуры в Maple. Для этого нужно просто переопределить параметры вывода командой interface( verboseproc=0,1,2,3)

• 0 - не выводить текст процедуры, только proc(x) ... end.
• 1 - выводить текст процедур пользователя (по умолчанию).
• 2 - выводить текст всех процедур.
• 3 - выводить текст и таблицу всех используемых процедур. 

 

Совет: если вы изменили значение переменной verboseproc, верните ей затем её старое значение 1.

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

 

Пример 6 (текст процедуры collect )

interface( verboseproc = 2,prettyprint=1,version );

print (collect);

interface( verboseproc = 1 );