Читать онлайн «Maple 9.5/10 в математике, физике и образовании». Страница 55

В пакете векторных операций определен ряд типовых операций с кривыми. Ниже представлено задание эллиптической кривой и вычисление в аналитической форме нормали и радиуса кривизны (файл vopcurves):

> SetCoordinates(cartesian);

cartesian

> assume(t::real);

> ell := <2*cos(t),sin(t)>;

ell := 2 cos(t)ex + sin(t)ey

> nv := simplify(PrincipalNormal(ell,t));

> len := simplify(LinearAlgebra:-Norm(nv, 2));

> r := simplify(RadiusOfCurvature(ell));

Теперь можно представить саму кривую (эллипс) и ее эволюту (рис. 4.39):

> ev := simplify(ell + r * nv / len);

> plot([[ell[1], ell[2], t=0..2*Pi], [ev[1], ev[2], t=0..2*Pi]]);

Рис. 4.39. Графики кривой — эллипса и ее эволюты

Нетрудно заметить, что для эллипса эволюта представляет собой удлиненную астроиду.

Для вычисления кривизны кривой С используется функция Curvature(C, t) в которой параметр t может и отсутствовать:

> Curvature(<cos(t),t,sin(t)>, t);

> с := Curvature(t -> <t,t^2,t^4>):

> simplify(c(t)) assuming t::real;

> SetCoordinates('polar');

polar

> Curvature(<exp(-t^2), t>):

> simplify(%) assuming t::real;

В аспекте практических приложений векторного анализа и теории поля особый интерес представляют приложения интегрирования пакете VectorCalculus. Так, видоизмененная функция int(f, dom) задает вычисление интеграла от функции f по области dom, например (файл vecint):

> restart:with(VectorCalculus):

> int(х^2+у^2, [x,y] = Circle(<0,1>, r));

> int(sin(х)*cos(у)*tan(z), [x,y,z] = Parallelepiped(0..Pi, 0..Pi/3, 0..Pi/4));

½√3 ln(2)

Функция PathInt(f, dom) вычисляет интеграл пути для функции f с Rn до R:

> PathInt(х^2, [х,y] = Line(<0,0>, <1,2>));

> PathInt(х^2+y^2, [х,y] = Circle(<0,0>, 3/2));

> PathInt(1, [х,y] = Ellipse(х^2+y^2/2-1));

Другая функция LineInt(F, dom), где F — вектор или процедура задания векторного поля, dom — параметр, характеризующий направление интегрирования, задает вычисление линейного интеграла в пространстве Rn:

> SetCoordinates(cartesian[х,y]);

cartesianx, у

> LineInt(VectorField(<х,y>), Line(<0,1>, <2,-5>));

14

> LineInt(VectorField(<y,-х>), Circlet<0,0>, r));

-2 r² π

> LineInt(VectorField(<y,-х>), Ellipse(х^2/4+y^2/9-1));

-12π

> LineInt(VectorField(<y,-х>), Arc(Ellipse(х^2/4+у^2/9-1), 0, Pi/2));

-3π

Функция ArcLength(C,dom) задает вычисление длины дуги С по известному интегральному выражению для нее:

> ArcLength(<r*cos(t),r*sin(t)>, t=0..Pi) assuming r>0;

πr

> ArcLength(t -> <t,t^2>, 0..2);

√17-¼ln(-4+√17)

> evalf(%);

4.646783762

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

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

Hessian(f, t) — создание матрицы гессиана;

Jacobian(f, v, det) — создание матрицы якобиана;

Wronskian(f, t) — создание матрицы вронскиана.

Примеры задания таких матриц приведены ниже (файл vecmatrix):

> Hessian(ехр(х*y), [х,y]);

> Hessian(а/(х^2+y^2+z^2), [х, y, z]);

> Н := unapply(%, [a,x,y,z]):

> Н(1/2, 0.3, 0.7, 0.1);

> Jacobian([r*cos(t), r*sin(t)], [r,t]);

> Jacobian([r*cos(t), r*sin(t)], [r,t], 'determinant');

> Wronskian([exp(t),ln(t),sin(t)], t);

> Wronskian([t, t^2, t^3], t)

К основным функциям теории поля относятся:

Curl(F) — вычисляет вихрь векторного поля в R³;

Divergence(F) — вычисляет дивергенцию векторного поля;

Flux(f, dom) — вычисляет поток векторного поля в R³;

Gradient(f, с) — вычисляет градиент функции f в пространстве от Rn до R;

Del(f, с) и Nabla(f, с) — векторные дифференциальные операторы;

Laplacian(f, с) или Laplacian(F) — вычисляет лапласиан функции f или векторного определения (процедуры) F;

ScalarPotential(v) — вычисляет скалярный потенциал векторного поля;

Torsion(C, t) — вычисляет торсион в R³;

VectorPotential(v) — вычисляет векторный потенциал в R³;

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

> restart:with(VectorCalculus): SetCoordinates('cartesian'[x,y,z]);

cartesianx, у, z

> F := VectorField( <-y,x,0> );

F:=-yēx +хēу

> Curl(F);

2ēz

> Del &x F;

2ēz

> Nabla &x F;

2ē

> CrossProduct(Del, F);

2ēz

> F := VectorField(<х^2,y^2,z^2>);

F:=-x²ēх + y²ēу + z²ēz

> Divergence(F);

2х + 2у + 2z

> Flux(VectorField(<x,y,z>, cartesian[x,y,z]), Sphere(<0,0,0>, r));

4r³ π

> Gradient(х^3/3+у^2, [x,y]);

x²ēx + 2yēу 0ēх

> Del(х^2+у^2+z^2);

2xēx + 2уēу + 2zēz

> Nabla(х^2+у^2+z^2);

2xēx + 2уēу + 2zēz

> Del . %;

6

> Laplacian(х^2+у^2+z^2, [x,y,z]);

6

> Laplacian(f(r,theta,z));

> SetCoordinates('cylindrical' [r, theta, z])

cylindricalr, θ, z

> Laplacian(f(r, theta, z));

> SetCoordinates('cartesian'[x,y,z]);

cartesianx, y, z

> v := VectorField(<x,y,-z>);

v := xēx + уēу - zēz

> ScalarPotential(v);

> v := VectorField(<-y,0,z>);

v := -yēx + zēz

> ScalarPotential(v); den := х^2 + y^2 + z^2;

den := x² + y² + z²

> ScalarPotential((x,y,z) -> <x,y,z>/den);

(x,y,z)→½ ln(x² + y² + z²)

> SetCoordinates('spherical'[r,phi,theta]);

sphericalr, φ, θ

> v := VectorField(<r,0,0>);

v:= r ēг

> ScalarPotential(v);

> restart:with(VectorCalculus): simplify( Torsion(<t,t^2,t^3>)) assuming t::real;

> Torsion(t -> <2*t,sin(t),cos(t)>);

> SetCoordinates('cartesian'[x,y,z]); v := VectorField(<y,-x,0>);

cartesianx, y, z v:= уēx - хēу

> VectorPotential(v);

-xzēx - yzēу

> SetCoordinates('cylindrical'[r,theta,z]);

cylindricalr, θ, z

> v := VectorField(<r,0,-2*z>);

v:= rēr -2zēz

> VectorPotential(v);

(-r sin(θ)² z - r cos(θ)² z) ēθ

> simplify(Curl(%));

rēr - 2zēz

Обратите внимание на то, что для гарантии правильного выполнения этих команд и отсутствия «зависания» компьютера может потребоваться команда restart и перезагрузка пакета VectorCalculus.

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

Рис. 4.40 Сложная поверхность с эффектами ее освещения внешним источником света

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

> J := Jacobian(f, [х, у, z]);

> J := DeleteColumn(J, [3]);

Тогда площадь поверхности вычисляется следующим образом:

> Int(Int(dA, x=0..2*Pi), y=0..2*Pi);

К сожалению, этот двойной интеграл Maple не вычисляет из-за сложности подынтегрального выражения, график которого представлен на рис. 4.41.

Рис. 4.41. График подынтегрального выражения

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