1610912777-ff63a1b83b9ac0b597c9346050946007 (824719), страница 8
Текст из файла (страница 8)
Примеры2172end;Процедура do_ident_mat(n) возвращает единичную матрицу размерности n.procedure do_ident_mat(n);beginscalar j1, j2;matrix mat_out(n,n);for j1:=1:n dofor j2:=1:n doif j1 = j2 thenmat_out(j1,j2):=1elsemat_out(j1,j2):=0;return mat_out;end;Процедура do_jordc_mat(n, lam) возвращает клетку Жордана размерности n и с собственным значением lam.procedure do_jordc_mat(n, lam);beginscalar j1;mat_out := do_ident_mat(n);mat_out := lam*mat_out;for j1 := 2:n domat_out(j1-1,j1) := 1;return mat_out;end;load_package normform;procedure solv_line_ode(equ, tt);beginscalar list_out,j1,j2,j3,j4,n,lam,u1,u2;n := part(equ,1);matrix amat(n,n), bmat(n,1), bmat2(n,1);for j1 := 1:n doГлава 15.
Примерыbeginu1 := part(equ,2*j1+1);for j2 := 1:n dobeginu2 := part(equ,2*j2);amat(j1,j2) := df(u1,u2);u1 := u1 - amat(j1,j2)*u2;end;bmat(j1,1) := u1;bmat2(j1,1) := mkid(part(equ,2*j1),1);end;s := jordan(amat);amat := part(s,1);s1 := get_jcell(s);bmat1 := part(s,3)*bmat;bmat2 := part(s,3)*bmat2;j4 := 0;for j1:=1:length(s1) dobeginj3 := part(s1,j1,1);lam := part(s1,j1,2);if j3 > 1 thenbeginmatrix amat3(j3,j3), bmat3(j3,1);amat3 := do_jordc_mat(j3, 0);amat4 := do_ident_mat(j3);amat5 := amat4;for j2 := 1:(j3 - 1) dobeginamat4 := amat4*amat3*tt/j2;amat5 := amat5 + amat4;end;for j2 := 1:j3 dobmat3(j2,1) := bmat2(j4 + j2,1);bmat3 := amat5*bmat3;if lam neq 0 thenbeginbmat3 := exp(repart(lam)*tt)*(cos(impart(lam)*tt)73Глава 15.
Примеры+ i*sin(impart(lam)*tt))*bmat3;end;matrix bmat4(j3,1);ind0 := 1;for j2 := 1:j3 dobeginbmat4(j3,1) := bmat1(j4 + j2,1);if bmat4(j3,1) neq 0 thenind0 := 0;end;if ind0 = 1 thenbeginamat3 := do_jordc_mat(j3, lam);if lam = 0 thenbeginbmat4 := bmat4*tt;bmat3 := bmat3 + bmat4;for j2 := 2:j3 dobeginbmat4 := amat3*bmat4*tt/j2;bmat3 := bmat3 + bmat4;end;endelsebmat3 := bmat3 - amat3^(-1)*bmat4;end;for j2 := 1:j3 dobmat2(j4 + j2,1) := bmat3(j2,1);endelsebeginj2 := j4 + 1;if lam = 0 thenbmat2(j2,1) := bmat2(j2,1) + bmat1(j2,1)*ttelsebmat2(j2,1) := exp(repart(lam)*tt)*(cos(impart(lam)*tt)+ i*sin(impart(lam)*tt))*bmat2(j2,1) + bmat1(j2,1)/lam;end;74Глава 15.
Примеры75j4 := j4 + j3;end;bmat2 := part(s,2)*bmat2;list_out := {};for j1 := 1:n dolist_out := append(list_out, {part(equ,2*j1) = bmat2(j1,1)});return list_out;end;15.5Точечные симметрии дифференциальныхуравненийВ данном разделе рассматриваются системы дифференциальных уравнений видаΦ(x, u, ∂x u, . . . , ∂xk u) = 0,x ∈ X = Rn , u : Rn → Y0 , Y0 = Rm ,∂xi u : X → Yi , Yi = Rm ⊗ S i Rn , i = 1, . .
. , k,Φ : Zk → Rs , Zi = X × Y0 × · · · × Yi , i = 0, 1 . . .(15.5.6)Групповой анализ дифференциальных уравнений использует геометрическое представление дифференциальных уравнений и их решений,как конечномерных многообразий в продолженном пространстве Zk .Здесь рассматриваются точечные преобразования, т.е. преобразования пространства зависимых X = Rn и независимых Y0 = Rm переменных, которые для нужд дифференциальных уравнений продолжаютсяна пространство полилинейных отображений до соответствующего порядка.Однопараметрическую локальную группу Ли точечных преобразований вполне определяет инфинитезимальный операторL=nXiξ (x, u)∂xi +i=1mXζ j (x, u)∂uj(15.5.7)j=1продолжение которого на пространство Zk записывается в виде|α|6kkL =L+X|α|>0ζα ∂uα ,(15.5.8)Глава 15.
Примеры76гдеζα+γj = Dj ζα − uα+γinXDj ξ i ,j = 1, . . . , n,(15.5.9)i=1Dj = ∂xj +Xuα+γj ∂uα ,j = 1, . . . , n.(15.5.10)|α|>0PЗдесь α = (α1 , . . . , αn ) — целочисленные мультииндексы, |α| = ni=1 αi , γj— мультииндекс, у которого j-я компонента равна единице, а остальныеравны нулю и ζ(0,...,0) = ζ.Система (15.5.6) допускает оператор (15.5.7), если и только если выполняется соотношениеLk ΦΦ=0 = 0.(15.5.11)Предлагаемый ниже алгоритм строит систему уравнений (15.5.11) ичасть простейших уравнений этой системы интегрирует.Функция compare_lists(s1, s2, n) возвращает 1, если совпадают первые n элементов списков s1 и s2.
В противном случае функция возвращает 0.procedure compare_lists(s1, s2, n);beginscalar z, j;z := 1;for j := 1:n dobeginif part(s1,j) neq part(s2,j) thenbeginz := 0;goto m1;end;end;m1:;return z;end;Функция copy_lists(s1, n) возвращает список составленный из первыхn элементов списка s1.Глава 15. Примеры77procedure copy_lists(s1, n);beginscalar s2, j;s2 := {};for j := 1:n dos2 := append(s2,{part(s1,j)});return s2;end;Функция do_ind(k, n) создает список всех мультииндексов длины nдо k-го порядка включительно.
Список содержит не только мультииндексы, но и информацию о том, какой операцией и из какого мультииндексапредыдущего порядка данный мультииндекс может быть получен. Длякаждого мультииндекса α = (α1 , . . . , αn ) порядка i существуют мультииндекс β порядка i − 1 и целое число j из интервала [1, n] такие, чтоα = β + γj (γj — мультииндекс, у которого j-я компонента равна единице, а остальные равны нулю). Мультииндекс β и число j определяютсянеоднозначно, но для реализации формулы продолжения (15.5.9) необходим какой-нибудь один вариант. Например, вызов функции do_ind(2, 2);возвратит список:{{{0,0,0,0}}, {{1,0,1,1},{0,1,1,2}}, {{2,0,1,1},{1,1,1,2},{0,2,2,2}}}.Здесь первых два элемента списков третьего уровня — мультииндексы, третий элемент — номер мультииндекса в списке мультииндексовпредыдущего порядка и четвертый элемент — число j.procedure do_ind(k, n);beginscalar s1, s2, s3, v, bb, j, k1, i1, i2;s1 := {};for j := 1:(n+2) dos1 := append(s1, {0});v := {{s1}};for k1 := 1:k dobegins1 := {};for j := 1:length(part(v,k1)) dobegins2 := copy_lists(part(v,k1,j), n);Глава 15.
Примеры78for i2 := 1:n dobegins3 := part(s2,i2) := part(s2,i2) + 1;bb := 1;for i1 := 1:length(s1) dobeginif compare_lists(part(s1,i1), s3, n) = 1 thenbeginbb := 0;goto m1;end;end;m1:;if bb = 1 thenbegins3 := append(s3,{j,i2});s1 := append(s1,{s3});end;end;end;v := append(v,{s1});end;return v;end;Функция do_var(u, indep, s) создает переменную из продолженногопространства Zk , используя зависимую переменную u, список независимых переменных indep и мультииндекс s.
Например, вызов функцииdo_var(u, {x, y}, {2, 1}); вернет переменную uxxy.procedure do_var(u, indep, s);beginscalar w, j, n, j1;n := length(indep);w := u;for j := 1:n dobeginfor j1 := 1:part(s,j) doГлава 15. Примеры79w := mkid(w,part(indep,j));end;return w;end;Функция prol_vars(indep, dep, l_ind) возвращает список переменныхпространства Zk , где indep — список независимых переменных, dep —список зависимых переменных и l_ind — список мультииндексов. Например, вызов функции prol_vars({x, y}, {u, v}, do_ind(2, 2)); вернет список:{x, y, u, v, ux, vx, uy, vy, uxx, vxx, uxy, vxy, uyy, vyy}.procedure prol_vars(indep, dep, l_ind);beginscalar s_out, n, m, s, s1, u, j, j1, j2, j3;s_out := indep;n := length(indep);m := length(dep);for j1 := 1:length(l_ind) dobegins := part(l_ind, j1);for j2 := 1:length(s) dobeginfor j := 1:m dobeginu := do_var(part(dep,j), indep, part(s, j2));s_out := append(s_out, {u});end;end;end;return s_out;end;Функция full_dif(indep, dep, l_ind) возвращает двухуровневый список, внутренние списки которого содержат компоненты векторов полногодифференцирования (15.5.10) (Dj , j = 1, .
. . , n), где indep — список независимых переменных, dep — список зависимых переменных и l_ind —Глава 15. Примеры80список мультииндексов. Например, команда full_dif({x, y}, {u, v}, do_ind(1, 2));вернет список:{{1, 0, ux, vx, uxx, vxx, uxy, vxy}, {0, 1, uy, vy, uxy, vxy, uyy, vyy}}.procedure full_dif(indep, dep, l_ind);beginscalar s_out, n, m, s, s1, s2, s3, u, j, j0, j1, j2;s_out := {};n := length(indep);m := length(dep);for j0 := 1:n dobegins3 := {};for j := 1:n doif j = j0 thens3 := append(s3, {1})elses3 := append(s3, {0});for j1 := 1:length(l_ind) dobegins := part(l_ind, j1);for j2 := 1:length(s) dobeginfor j := 1:m dobegins1 := copy_lists(part(s,j2), n);s2 := part(s1,j0) := part(s1,j0) + 1;u := do_var(part(dep,j), indep, s2);s3 := append(s3, {u});end;end;end;s_out := append(s_out, {s3});end;return s_out;end;Глава 15.
Примеры81Функция field_act(s1, s2, u) возвращает результат действия векторного поля на выражение u, где s1 — список переменных и s2 — списоккомпонент поля.procedure field_act(s1, s2, u);beginscalar s_out, j;s_out := 0;for j := 1:length(s1) dobeginif part(s2, j) neq 0 thenbegins_out := s_out + part(s2,j)*df(u, part(s1,j));end;end;return s_out;end;Функция prol_al(s_in, indep, dep, l_ind, vars, full_d), применяя формулы (15.5.9), возвращает список компонент продолженного векторного поля (15.5.8), где s_in — компоненты исходного векторного поля L,indep — список независимых переменных, dep — список зависимых переменных, l_ind — список мультииндексов, vars — список переменныхпродолженного пространства Zk , full_d — список компонент операторовполного дифференцирования.procedure prol_al(s_in, indep, dep, l_ind, vars, full_d);beginscalar s_out, n, m, m1, m2, m21, m3, m4, s, s0, s2, s3,j, j1, j2, j3, u, u1, u2;s_out := s_in;n := length(indep);m := length(dep);m1 := n;m21 := m;for j1 := 2:length(l_ind) dobeginm2 := 0;Глава 15.
Примеры82s := part(l_ind, j1);for j2 := 1:length(s) dobeginm3 := part(s,j2,n+1);m4 := part(s,j2,n+2);for j := 1:m dobeginu1 := part(s_out, m1 + (m3-1)*m + j);s0 := part(full_d,m4);u := field_act(vars, s0, u1);s2 := copy_lists(part(l_ind,j1-1,m3), n);for j3 := 1:n dobegins3 := part(s2,j3) := part(s2,j3) + 1;u2 := do_var(part(dep,j), indep, s3);u := u - u2*field_act(vars, s0, part(s_in, j3));end;m2 := m2 + 1;s_out := append(s_out, {u});end;end;m1 := m1 + m21;m21 := m2;end;return s_out;end;Процедура all_coeff(ur, var, j) находит все коэффициенты полиномовиз списка ur от переменных из списка var и выводит их на консоль илив файл.
Параметр j играет роль счетчика переменных при рекурсивном вызове процедуры. Перед вызовом функции all_coeff должна бытьинициализирована глобальная переменная neqw, в которой сохраняетсяобщее количество коэффициентов.procedure all_coeff(ur, var, j);beginscalar ur1, v, i1, i2, i3;for i1 := 1:length(ur) doГлава 15. Примеры83beginur1 := coeff(part(ur, i1), part(var, j));if j < length(var) thenall_coeff(ur1, var, j + 1)elsebeginfor i2:=1:length(ur1) dobeginv := part(ur1, i2);if v neq 0 thenbeginneqw := neqw + 1;write "eqw(", neqw, "):=", v,"$";end;end;end;end;end;Функция no_dep(eqw, nn) из массива дифференциальных уравненийeqw выбирает уравнения первого порядка, в которых всего одно слагаемое. Такие уравнения означают, что определенная компонента поля не зависит от определенной переменной (например, уравнение −3 ∗df (ky, x) = 0) означает, что компонента поля ky не зависит от переменнойx. Функция no_dep возвращает двухуровневый список.
Двухэлементныесписки второго уровня содержат компоненты поля и переменные, от которых эти компоненты не зависят.procedure no_dep(eqw, nn);beginscalar s_out, ind, bb, ob, i, j, k, z;bb := 1;ob := 0;s_out := {};while bb eq 1 dobeginob := ob + 1;bb := 0;Глава 15. Примеры84for i := 1:nn dobeginif eqw(i) neq 0 thenbeginif part(eqw(i),0) eq minus theneqw(i):=-eqw(i);ind := 0;if (part(eqw(i),0) eq df) and (arglength(eqw(i)) < 3) thenbeginind := 1;z := eqw(i);end elseif part(eqw(i),0) eq times thenbeginfor k := 1:arglength(eqw(i)) dobeginif (arglength(part(eqw(i), k)) > 0) and(part(part(eqw(i), k),0) eq df) and(arglength(part(eqw(i), k)) < 3) thenbeginz := part(eqw(i), k);ind := ind + 1;end;end;end;if ind eq 1 thenbeginbb := 1;s_out := append(s_out,{{part(z,1),part(z,2)}});for j := 1:nn dobegineqw(j) := sub({z=0}, eqw(j));end;end;end;end;end;return s_out;Глава 15.