Мостовской А.П. Информационные технологии в математике

Подождите немного. Документ загружается.

1.19 Интерполяция данных 61

Trapezoidal - назначает для решения рекурсивный метод трапеций.

Он особенно успешен, если подынтегральная функция периодична и

интервал интегрирования составляет точно один период. Для много-

мерного интегрирования данный метод обращается к декартову произ-

ведению одномерных правил трапеции.

MultiDimensional - назначает для исполнения многомерный алго-

ритм, что имеет смысл только для кратных интегралов.

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

Гаусса-Кронрода. Дополнительная опция GaussPoints−> n

устана-

вливающая количество точек в Гауссовой части квадратуры Гаусса-

Кронрода. По умолчанию значение этой опции есть Automatic.

SingularityDepth−> n

опция указывает число рекурсий, приме-

няемых перед тем, как произойдет изменение переменной на границах

интервала интегрирования. По умолчанию n

= 4.

MinRecursion−> 0, MaxRecursion−> 6 - минимальное и макси-

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

Приведены значения по умолчанию.

Примеры.

1). NIntegrate[Sqrt[2*x+1],x,0,1] 1.39872

2). NIntegrate[1/(x*y),x,4,4.4,y,2,2.6] 0.025006

3). NIntegrate[x ∗ y, x, 0, 1, y, x, x

, z, x ∗ y, x

∗ y

] 0.010582

4). NIntegrate[1/Sqrt[1 − x

], x, 0, 1] 1.21433

5). NIntegrate[Exp[−x

], {x, −1000, 1000}, MinRecursion−> 3,

MaxRecursion−> 10] ⇒ 1.77245

6). NIntegrate[Exp[−x

], {x, −1000, −10, 10, 1000}] ⇒ 1.77245

1.19 Интерполяция данных

Решении некоторых задач ( например, численное решение диффе-

ренциальных уравнений ) получается в виде ряда точек. Математика

позволяет находить интерполяционные многочлены таких решений.

Команда

Interpolation[data], (1)

где data = {{x

, f

}, {x

, f

}, ...{x

, f

}}, создает функцию, которая при-

нимает в точках x

значение f

. Данные data можно представить в

форме {f

, f

, ..., f

}, при этом считается, что x

= i, i = 1, 2, ...n.

62 ГЛАВА 1. СИСТЕМА МАТЕМАТИКИ

Выполнение команды (1) завершается выводом вида

Out[m] = InterpolatingFunction[{{x

, x

}}, <>].

Для дальнейшей работы с этой функцией надо обозначить f=%m.

Можно команду (1) написать в виде f= Interpolation[data] и дальше

производить вычисления, например, найти f[2] и так далее.

Пример. Пусть

p = T able[{x, 1/x}, {x, 1, 10}] ⇒ {{1, 1}, {2, 1/2}, ..., {10, 1/10}} - спи-

сок из координат десяти точек гиперболы. Построим интерполяцион-

ный многочлен:

f = Interpolation[p] ⇒ InterpolatingF unction[{{1, 10}}, <>].

Теперь можно найти, например, f[3.] ⇒ 0.333333 или построить график

разности P lot[1/x − f[x], {x, 1, 10}, P lotRange−> All].

Команда

InterpolatingPolynomial[data, var]

возвращает интерполяционный многочлен от переменной var для точек

с координатами, принадлежащими списку

data = {{x

, y

}, {x

, y

}, ..., {x

, y

}}.

Многочлен принимает значения y

при var = x

. Список можно пред-

ставить в виде data = {y

, y

, ..., y

}, при этом координаты x точек

предполагаются равными 1,2,...,n соответственно.

Пример.

InterpolatingP olynomial[{1, 3, 6}, x] ⇒ 1 + (2 +

(−2 + x)(−1 + x).

%//Simplify ⇒

x(x+1)

- парабола, проходящая через точки (1,1),

(2,3), (3,6).

1.20 Решение дифференциальных уравнений

Математика имеет возможность решать обыкновенные линейные

дифференциальные уравнения 1-го и 2-го порядка, некоторые нели-

нейные дифференциальные уравнения, системы линейных дифферен-

циальных уравнений, уравнения и системы уравнений в частных про-

изводных

Для аналитического решения дифф еренциального уравнения приме-

няется команда

DSolve[eqns, y[x], x],

Последнее относится к возможностям Mathematica 5.

1.20 Решение дифференциальных уравнений 63

список eqns содержит дифференциальное уравнение относительно функ-

ции y[x] и, при необходимости, начальные условия, записанные как

уравнения: y[x

] == y

, y

] == y

, .... Если eqns состоит только из

дифференциального уравнения, то фигурные скобки, обозначающие спи-

сок, можно не писать. Решение дифференциального уравнения записы-

вается в виде {{y[x]−> expr}}, означающем, что решением являет-

ся функция y[x] = expr. При этом функция expr может быть чистой

функцией одной переменной или интерполяционным многочленом.

Команда

DSolve[{eqn

, eqn

, ..., eqn

, ..., x

], ..., y

, ..., x

]}, {x

, ..., x

}]

решает систему дифференциальных уравнений. В список могут вхо-

дить и начальные условия.

Опция DSolveConstants−> K в DSolve заменяет постоянную ин-

тегрирования, традиционно обозначаемую буквой C, на букву K.

Примеры.

1). DSolve[y

[x] == 2 a x

, y[x], x] ⇒ {{y[x]−>

+ C[1]}}.

2). DSolve[{y[x] == z

[x], z[x] == −y

[x]}, {y[x], z[x]}, x] ⇒

{{y[x]−> C[1]Cos[x] − C[2]Sin[x], z[x]−> C[2]Cos[x] + C[1]Sin[x]}}.

3). s = DSolve[{y

[x] + y[x] == x, y[0] == 1}, y[x], x] ⇒ {{y[x]−>

−x

(2 −E

+ E

x}}. Построим график полученного частного решения:

Plot[y[x]/.%,{x,-10,10}].

Можно построить этот график, не используя решение уравнения s:

Plot[Evaluate[y[x]/.s],{x,-10,10}].

4). DSolve[y

[x] − y

[x] − 6 y[x] == 0, y[x], x, DSolveConstants− >

K] ⇒ {{y[x]−> E

−2x

K[1] + E

K[2]}}.

Для численного решения дифференциальных уравнений использует-

ся команда:

NDSolve[eqns, y, {x, x

min

, x

max

}],

где список eqns содержит дифференциальное уравнение относительно

неизвестной функции y[x], x ∈ [x

min

, x

max

] и начальные условия, за-

писанные в виде y[x

] == y

, .... Решение такого дифференциального

уравнения получается в виде интерполяционной функции.

Команда

NDSolve[eqns, {y

, y

, ...}, {x, x

min

, x

max

}]

ищет численные решения системы дифференциальных уравнений от-

носительно функций y

[x] при x ∈ [x

min

, x

max

MaxSteps−> m опция к NDSolve, которая определяет максималь-

64 ГЛАВА 1. СИСТЕМА МАТЕМАТИКИ

ное количество шагов m.

StartingStepSize−> w опция к NDSolve, определяющая величину

w начального шага.

Примеры.

1). NDSolve[{y

[x] == y[x], y[0] == 1}, y, {x, 0, 2}] ⇒

{{y−> InterpolatingF unction[{{0., 2.}}, <>]}}.

Можно построить график решения: Plot[y[x]/.%,{x,0,2}].

2). Пусть

dif = {x

[t] == −y[t] − x[t]

, y

[t] == 2 x[t] − y[t], x[0] == y[0]] == 1}.

Решаем систему:

sol = NDSolve[dif, {x, y}, {t, 0, 10}] ⇒

{{x−> InterpolatingF unction[{{0., 10.}}, <>],

y−> InterpolatingF unction[{{0., 10.}}, <>]}}.

Найдем значения решений при t=3:

{x[3], y[3]}/.sol ⇒ {{−0.13973, −0.517751}}.

Построим график решения:

ParametricPlot[Evaluate[{x[t],y[t]}/.sol],{t,0,10}].

3). Орбита планет. Траектория точечной массы в гравитационном

поле большей массы , помещенной в начало координат, описывается си-

стемой дифференциальных уравнений

= −

где r =

+ y

Найдем численное решение этой системы:

diff = NDSolve[

[t] == −x[t]/(x[t]

+ y[t]

)

[t] == −y[t]/(x[t]

+ y[t]

)

x[0] == 1, x

[0] == 0.2, y[0] == 0, y

[0] == 1.25},

{x, y}, {t, 0, 45}] ⇒

{{x−> InterpolatingF unction[{{0., 45.}}, <>],

y−> InterpolatingF unction[{{0., 45.}}, <>]}}

Построим решение:

P arametricP lot[Evaluate[{x[t], y[t]}/.diff], {t, 0, 45}].

4). Решим систему трех дифференциальных уравнений:

NDSolve[{x

[t] == −3(x[t] − y[t]), y

[t] == −x[t] z[t] + 27 x[t] − y[t],

[t] == x[t] y[t] − z[t], x[0] == z[0] == 0, y[0] == 1},

{x, y, z}, {t, 0, 20}, MaxSteps−> 3000] ⇒

1.21 Преобразование Лапласа 65

{{x−> InterpolatingF unction[{{0., 20.}}, <>],

y−> InterpolatingF unction[{{0., 20.}}, <>],

z−> InterpolatingF unction[{{0., 20.}}, <>], }}.

Построим график решения:

P arametricP lot3D[Evaluate[{x[t], y[t], z[t]}/.%], {t, 0, 20},

P lotP oints−> 1000]

Задачи. Решить дифференциальные уравнения.

1. y

= y tg(x).

2. y

− y tn(x) = sec(x).

3. y

− 2xy = 1.

4. x

+ 3xy = (sin x)/x.

5. y

+ x y

+ y = 0.

Задачи. Решить задачи Коши. Построить графики решений.

1. −x

+ y

+ 3xy + x

= 0; y(1) = 0.

2. (x

+ sin x + 2)y

+ 2xE

+ y cos x = 0; y(0) = 1.

3. x

+ xy(x y + 4) + 2 = 0; y(2) = 0.

4. xy

+ axy

+ 2y + bx = 0; y(1) = 0, при a = 1, b = 4.

5. y

= y

; y(1) = 1, y

(1) = −1.

1.21 Преобразование Лапласа

Команды, связанные с преобразованием Лапласа, становятся доступ-

ными после загрузки пакета расширения

<< Calculus`LaplaceTransform`

66 ГЛАВА 1. СИСТЕМА МАТЕМАТИКИ

Команда

LaplaceTransform[f[t], t, p, optins−> value]

сформирует функцию аргумента p, которая является преобразованием

Лапласа функции f[t]:

LaplaceTransform[f[t], t, p] =

Integrate[Exp[−p t]f[t], {t, 0, Inﬁnity}].

Стандартные примеры:

{LaplaceTransform[1, t, p],

LaplaceTransform[t, t, p],

LaplaceTransform[E

(a t)

, t, p],

LaplaceTransform[t

, t, p],

LaplaceTransform[Cos[w t], t, p],

LaplaceTransform[Cosh[w t], t, p], } ⇒

{

, p

−2

−a + p

, p

−1−n

Gamma[1 + n],

+ w

− w

Команда

InverseLaplaceTransform[g[p], p, t, optins−> value]

формирует функцию аргумента t, преобразованием Лапласа которой

является функция g[p].

Пример 1. Решим линейное дифференциальное уравнение с посто-

янными коэффициентами, правая часть которого есть функция , пре-

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

eq = y

[t] − 3y

[t] + 2y[t] == 4t + E

Преобразование Лапласа преобразует это дифференциальное уравне-

ние в алгебраическое относительно преобразования Лапласа функции

y[t]:

Leq = LaplaceT ransform[eq, t, p] ⇒

2 LaplaceT ransform[y[t], t, p] + s

LaplaceT ransform[y[t], t, p]−

−3(p LaplaceT ransform[y[t], t, p] −y[0]) −p y[0] −y

[0] ==

−3 + p

Решаем это алгебраическое уравнение относительно

LaplaceT ransform[y[t], t, p] :

Aeq = Solve[Leq, LaplaceT ransform[y[t], t, p]] ⇒

{{LaplaceT ransform[y[t], t, p] → −((12−4p−p

−9p

y[0]+6p

y[0]−p

y[0]+

1.21 Преобразование Лапласа 67

[0] − p

[0])/(−6p

+ 11p

− 6p

+ p

))}}.

Применяя теперь команду InverseLaplaceTransform, получаем решение

дифференциального уравнения eq:

Ans = InverseLaplaceT ransform[Aeq[[1, 1, 2]], p, t] ⇒

3 +

+ 2t +

(−7 + 4y[0] − 2y

[0])

+ E

(−y[0] + y

[0]).

Сделаем проверку - подставим найденную функцию и ее производные

в уравнение eq:

eq/.{y[t]−> Ans, y

[t]−> D[Ans, t], y

[t]−> D[Ans, {t, 2}]}//Simplify ⇒

T rue.

Построим график решения дифференциального уравнения при на-

чальных условиях y[0]=1, y’[0]=1:

P lot[Ans/.{y[0]−> 1, y

[0]−> 1}, {t, 0, 1}].

Задачи. Применяя преобразование Лапласа, решить следующие диф-

ференциальные уравнения и построить графики решений при различ-

ных начальных условиях.

1. y

(t) − 2y(t) = sin t.

2. y

(t) − 4y

(t) + 4y(t) = t

3. y

(t) − 5y

(t) + 4y(t) = E

4. y

(t) + 2y

(t) + 2y(t) = t.

5. y

(t) + 2y

(t) + 2y(t) = E

−t

sin t.

6. y

000

(t) + y(t) = 0.

Пример 2. Рассмотрим систему двух линейных дифференциальных

уравнений с постоянными коэффициентами:

eq = {x

[t] == 2(y[t] − 2x[t]), y

[t] == 2(x[t] − 2y[t])}

и начальными условиями:

nu = {x[0]−> 1, y[0]−> 1, x

[0]−> Sqrt[3]Sqrt[2],

[0]−> −Sqrt[3]Sqrt[2]};

Решим систему eq. Для этого применим преобразование Лапласа к

системе eq:

lteq = LaplaceTransform[eq, t, p];

Разрешим систему алгебраических уравнений lteq относительно не-

известных LaplaceTransform[x[t], t, p], LaplaceTransform[y[t], t, p]. При

этом учитываются начальные условия:

aleq = Solve[lteq, {LaplaceT ransform[x[t], t, p],

68 ГЛАВА 1. СИСТЕМА МАТЕМАТИКИ

LaplaceT ransform[y[t], t, p]}]/.nu//Simplify;

Выделим функции - решения системы:

x = InverseLaplaceT ransform[aleq[[1, 1, 2]], p, t];

y = InverseLaplaceT ransform[aleq[[1, 2, 2]], p, t];

(1)

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

тыми в концах команд.

Сделаем проверку:

eq/.{x[t]−> x, y[t]−> y, x

[t]−> D[x, {t, 2}],

[t]−> D[y, {t, 2}]}//Simplify ⇒ {T rue, T rue}

Построим решения, то есть построим кривую, заданную параметри-

ческими уравнениями (1):

P arametricP lot[Evaluate[{x[t], y[t]}/.

{x[t]−> x, y[t]−> y}], {t, 0, 20}]

Задача. Применяя преобразование Лапласа, решения систему диф-

ференциальных уравнений с начальными условиями x(0)=2, y(0)= 1,

x’(0)= 3, y’(0)=1. Построить график решения.

(

(t) = −3x(t) + 4y(t) + cos t,

(t) = −2x

(t) + 3y(t) + t.

1.22 Функции минимизации и максимизации

Математика дает разнообразные возможности по решению задач

оптимизации - от поиска элементов списка с минимальным или мак-

симальным значением, до поиска локальных и даже глобальных мини-

мумов функций, заданных аналитически.

Для поиска максимального и минимального значений ряда величин,

входящих в список, применяются следующие команды:

Max[x

, x

, ...]

выбирает наибольшее из xi,

Max[{x

, x

, ...}, {y

, ...}, ...]

находит наибольший элемент среди элементов всех списков,

Min[x

, x

, ...]

находит наименьшее из x

Min[{x

, x

, ...}, {y

, ...}, ...]

выбирает наименьший элемент среди элементов всех списков.

Примеры.

1.22 Функции минимизации и максимизации 69

Max[1, 5, 2, 6.5, 3, 4] ⇒ 6.5, Max[{1, 3, 2}, {4, 5, 6}, {9, 8, 7}] ⇒ 9,

Min[1, 5, 2, 6.5, −3, 4] ⇒ −3, Min[{1, 3, 2}, {456}, {9, 8, 7}] ⇒ 1.

Для нахождения локального минимума функции применяется коман-

да

FindMinimum[f[x], {x, x

}],

- поиск локального минимума функции f[x], начиная со значения x = x

Ответ записывается в виде списка {f

, {x−> x

}}, где x = x

точка,

в которой функции f[x] достигает своего минимального значения f

Команда

FindMinimum[f[x], {x, x

, a, b}]

находит локальный минимум функции f[x] на интервале [a,b]. Началь-

ное значение x

∈ [a, b].

Команда

FindMinimum[f[x

, x

, ..., x

], {x

, a

}, {x

, a

}, ..., {x

, a

}]

находит локальный минимум функции f[x

, x

, ..., x

]. Начальное зна-

чение определяется точкой {a

, a

, ..., a

}. Ответ записывается в виде

списка

, {x

−> b

, x

−> b

, ..., x

−> b

где f

- минимальное значение функции, {b

, b

, ..., b

} - точка, в ко-

торой функция достигает данного минимума .

Примеры.

F indMinimum[−x Exp[−2 x], {x, 1}] ⇒ {−0.18394, {x−> 0.5}}.

F indMinimum[−x Exp[−2 x], {x, 0.2, 0, 1}] ⇒ {−0.18394, {x−> 0.5}}.

F indMinimum[−5 x Exp[−x/2] (2 + Sin[3 x]), x, 1] ⇒

{−7.17833, {x−> 0.783139}}.

F indMinimum[100 (y − x

)

+ (1 − x)

, {x, 0}, {y, 0}] ⇒

{9.98756 · 10

−13

, {x−> 1., y−> 0.999999}}.

Следующие две команды применяются для поиска глобального мак-

симума и минимума линейной функции в Mathematica 4 или в более

ранней версии:

ConstrainedMax[f[x

, x

, ..., x

], {ineqs}, {x

, x

, ..., x

}]

находит глобальный максимум, а

ConstrainedMin[f[x

, x

, ..., x

], {ineqs}, {x

, x

, ..., x

}],

находит глобальный минимум функции f[x

, x

, ..., x

] в области, опре-

деляемой неравенствами ineqs. Предполагается, что все переменные

, x

, ..., x

неотрицательны.

Пример. ConstrainedMin[x + 3y + 7z, {x − 3y < 7, 2x + 3z >= 5,

x + y + z < 10}, {x, y, z}] ⇒ {

, {x−>

, y−> 0, z−> 0}}.

70 ГЛАВА 1. СИСТЕМА МАТЕМАТИКИ

Постройте многоугольник-решение системы данных неравенств.

В Mathematica 5 для нахождения глобальных точек максимума и

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

Minimize, Maximize. Так, команда

Minimize[f[x, y, ...], {x, y, ...}]

находит глобальный максимум не обязательно линейной функции от-

носительно ее переменных, а команда

Minimize[f[x, y, ...], {conds}, {x, y, ...}]

решает аналогичную задачу, при условии, что переменные удовлетво-

ряют условиям conds.

Например,

Minimize[x

+3y +7z

, {x−3y < 7, 2x+3z >= 5, x+ y +z < 10, x >=

0, y >= 0, z >= 0}, {x, y, z}] ⇒ {

175

, {z− >

, x− >

, y− > 0}}.

В ответе указано минимальное значение функции и координаты точки,

в которой минимум достигается.

Задача. Для изготовления изделий a и b склад может отпустить

металла не более 80 кг, причем на изделие a расходуется 2 кг, а на b -

1 кг металла. Спланировать производство так, чтобы была обеспечена

наибольшая прибыль, если изделий a требуется изготовить не более 30

штук по цене 5 рублей за штука, а изделий b не более 40 штук по цене

3 рубля за штука.

1.23 Специальные математические функции

В ядро Математики включен ряд специальных математических функ-

ций. Приведем некоторые из них.

ChebyshevT[n, x], ChebyshevU[n, x], HermiteH[n, x], JacobiP[n,

a,b, x], GegenbauerC[n, m, x], LaguerreL[n, x] -

ортогональные многочлены Чебышева n-ой степени I рода, Чебыше-

ва n-ой степени II рода, Эрмита n-ной степени, Якоби n-ной степени,

Гегенбауэра, Лагерра n-ной степени соответственно.

AiryAi[z], AiryBi[z] - функции Эйри Ai(z) и Bi(z).

BesselJ[n, z], BesselY[n, z] - функция Бесселя I и II рода.

BernoulliB[n] - n-ое число Бернулли.

BernoulliB[n, x] - полином Бернулли n-ой степени.

Beta[a, b] - Эйлерова бета-функция B(a, b).

Binomial[n, m] - биномиальный коэффициент.

EllipticE[m] - полный эллиптический интеграл E(m).