Лабораторный практикум - Математические методы и модели в расчётах на ЭВМ
Подождите немного. Документ загружается.
Для схемы рис. 2 (φ
3
= 0; дерево графа ветви 1 и 5; ветви связи – 2, 3, 4;
главные контуры – ветви 2 и 1; 3, 1 и 5; 4 и 5) матрицы, входящие в эти
уравнения, имеют вид
;
11000
10101
00011
;
11100
00111
BA
)./1/1/1/1/1();(
;000;0
5432154321
435321
RRRRRdiagGRRRRRdiagR
JJJEEEEE
TT
Пусть Е
1
= 3 В, Е
2
= Е
3
= 2 В, Е
5
= 1 В; J
3
= 1 A, J
4
= 0,5 A;
R
1
= 4 Ом, R
2
= 6 Ом, R
3
= 2 Ом, R
4
= 8 Ом, R
5
= 10 Ом.
Рассмотрим решение матричных уравнений (1) и (2) в среде MATLAB.
>> A=[-1 1 -1 0 0;0 0 1 -1 1]
A =
-1 1 -1 0 0 - матрица соединений
0 0 1 -1 1
>> B=[1 1 0 0 0;-1 0 1 0 -1;0 0 0 1 1]
B =
1 1 0 0 0
-1 0 1 0 -1 - матрица главных контуров
0 0 0 1 1
>> RB=[4 6 2 8 10]
RB =
4 6 2 8 10 - матрица сопротивлений ветвей
>> R=diag(RB)
R =
4 0 0 0 0
0 6 0 0 0 - диагональная матрица сопротивлений
0 0 2 0 0
0 0 0 8 0
0 0 0 0 10
>> RK=B*R*B'
RK =
10 -4 0
-4 16 -10 - матрица контурных сопротивлений
0 -10 18
>> E=[3;2;2;0;1]
E =
3
2 - матрица ЭДС ветвей
2
0
1
31
>> J=[0;0;-1;-0.5;0]
J =
0
0 - матрица источников токов
-1.0000
-0.5000
0
>> EK=B*E-B*R*J
EK =
5 - матрица контурных ЭДС
0
5
>> IB=RK^(-1)*EK - токи ветвей связи
IB =
0.7161 ток I
2
0.5402 ток I
3
0.5779 ток I
4
>> IB=RK\EK - токи ветвей связи
IB =
0.7161
0.5402
0.5779
>> IB=inv(RK)*EK - токи ветвей связи
IB =
0.7161
0.5402
0.5779
>> I=B'*IB - токи обобщённых ветвей
I =
0.1759 ток I
1
0.7161
0.5402
0.5779
0.0377 ток I
5
Токи сопротивлений рассчитываются по формулам:
I
R1
= I
2
– I
3
= 0,716 – 0,540 = 0,156 A;
I
R2
= I
2
= 0,716 A;
I
R3
= I
3
– J
3
= 0,540 – 1 = - 0,46 A;
I
R4
= I
4
– J
4
= 0,578 - 0,5 = 0,078 A;
I
5
= I
4
– I
3
= 0,578 -0,540 = 0,038 A.
32
Расчёт потенциалов узлов:
>> GB=[1/4 1/6 1/2 1/8 1/10]
GB =
0.2500 0.1667 0.5000 0.1250 0.1000 - проводимости ветвей
>> G=diag(GB)
G =
0.2500 0 0 0 0
0 0.1667 0 0 0 - диагональная матрица
0 0 0.5000 0 0 проводимостей ветвей
0 0 0 0.1250 0
0 0 0 0 0.1000
>> GY=A*G*A'
GY =
0.9167 -0.5000 - матрица узловых
-0.5000 0.7250 проводимостей
>> JY=A*J-A*G*E
JY =
2.4167 - матрица узловых токов
-1.6000
>> fu=GY^(-1)*JY - матрица узловых потенциалов
fu =
2.2965 потенциал первого узла
-0.6231 потенциал второго узла
>> fu=GY\JY
fu =
2.2965
-0.6231
>> fu=inv(GY)*JY
fu =
2.2965
-0.6231
Таким образом, если А – матрица коэффициентов системы уравнений
размерностью n×n, В – вектор – столбец свободных членов размерностью
n×1, х – вектор – столбец неизвестных размерностью n×1, то решение
системы уравнений А*х = В возможно по одному любому из следующих
выражений: x = inv(A)*B; x =A\B и x = A
-1
*B при условии, что det(A) ≠ 0.
Решение систем линейных уравнений возможно также с помощью
функций solve ( ). Например, требуется решить систему уравнений
.5383
;462
;935
321
321
321
xxx
xxx
xxx
33
>> syms x1 x2 x3; - определение символьных переменных
>> X=solve('5*x1+3*x2-x3=9','2*x1-x2+6*x3=4','3*x1+8*x2+3*x3=5')
X =
x1: [1x1 sym]
x2: [1x1 sym]
x3: [1x1 sym]
>> X.x1 - вывод решения
ans =
13/7
>> X.x2
ans =
-10/119
>> X.x3
ans =
4/119
>> vpa(X.x1,4) - вывод решения с заданным числом знаков ответа
ans =
1.857
>> vpa(X.x2,4)
ans =
-.8403e-1
>> vpa(X.x3,4)
ans =
.3361e-1
Если матрица А – разрежённая матрица, то целесообразно применять
итерационные методы (метод наименьших квадратов, сопряжённых
градиентов и др.) решения таких систем уравнений и соответствующих им
функций. Эти же функции применяются и для решения обычных систем
линейных алгебраических уравнений:
bccg( ) – двунаправленный метод сопряжённых градиентов
>> x1=bicg(A,B); - выдаёт результат с погрешностью по умолчанию 10
-6
bicg converged at iteration 3 to a solution with relative residual 8.1e-016
>> x2=bicg(A,B,0.0001); - выдаёт результат с заданной погрешностью
bicg converged at iteration 3 to a solution with relative residual 8.1e-016
>> x3=bicg(A,B,0.0001,15); - выдаёт результат с заданной погрешностью
при заданном числе итераций
bicg converged at iteration 3 to a solution with relative residual 8.1e-016
>> X=[x1,x2,x3]
X =
1.8571 1.8571 1.8571
-0.0840 -0.0840 -0.0840
0.0336 0.0336 0.0336
во всех случаях получен одинаковый ответ;
cgs( ) – квадратичный метод сопряжённых градиентов
34
>> cgs(A,B)
cgs converged at iteration 3 to a solution with relative residual 3.3e-015
ans =
1.8571
-0.0840
0.0336
pcg( ) – метод сопряжённых градиентов
>> pcg(A,B)
pcg stopped at iteration 3 without converging to the desired tolerance 1e-006
because the maximum number of iterations was reached.
The iterate returned (number 2) has relative residual 0.44
ans =
2.0856
0.3726
-0.0385
bicgstab( ) – устойчивый двунаправленный метод
>> bicgstab(A,B)
bicgstab converged at iteration 2.5 to a solution with relative residual 1.8e-016
ans =
1.8571
-0.0840
0.0336
lsqr( ) – метод наименьших квадратов
>> lsqr(A,B)
lsqr stopped at iteration 3 without converging to the desired tolerance 1e-006
because the maximum number of iterations was reached.
The iterate returned (number 3) has relative residual 1.6e-015
ans =
1.8571
-0.0840
0.0336
gmres( ) – метод минимизации обобщённой невязки
>> gmres(A,B)
gmres converged at iteration 3 to a solution with relative residual 8e-017
ans =
1.8571
-0.0840
0.0336
35
qmr( ) – квазиминимизация невязки
>> qmr(A,B)
qmr converged at iteration 3 to a solution with relative residual 2e-016
ans =
1.8571
-0.0840
0.0336
Более подробное описание этих функций приведено в [4].
Эти же функции можно применять при расчётах электрических цепей
символическим методом (решение систем линейных алгебраических
уравнений с комплексными коэффициентами):
>> A=[1+2i -4+4i 0+i;-4-10i 5+5i 10+i;5+i 7+5i 9-5i]
A =
1.0000 + 2.0000i -4.0000 + 4.0000i 0 + 1.0000i
-4.0000 -10.0000i 5.0000 + 5.0000i 10.0000 + 1.0000i
5.0000 + 1.0000i 7.0000 + 5.0000i 9.0000 - 5.0000i
>> B= [100+60i;220-100i;-50i]
B =
1.0e+002 *
1.0000 + 0.6000i
2.2000 - 1.0000i
0 - 0.5000i
>> X=A\B
X =
-7.1841 +17.7280i
-11.2561 -26.7692i
-5.4121 + 9.4609i
>> lsqr(A,B)
lsqr stopped at iteration 3 without converging to the desired tolerance 1e-006
because the maximum number of iterations was reached.
The iterate returned (number 3) has relative residual 5.1e-015
ans =
-7.1841 +17.7280i
-11.2561 -26.7692i
-5.4121 + 9.4609i
Порядок выполнения работы
1. Для электрической цепи, заданной в виде графа (рис. 5), записать матрицы
соединений, главных контуров, источников ЭДС и токов, диагональные
матрицы сопротивлений и проводимостей ветвей. Номиналы элементов схем
приведены в таблице 1, 2.
2. Решить матричные уравнения метода контурных токов и узловых
потенциалов с помощью различных операторов. Токи сопротивлений
проверить по уравнению энергетического баланса.
36
Вопросы для самопроверки
1. Приведите примеры линейных элементов электрических цепей и их
математических моделей.
2. Как формируются матрицы инциденций и главных контуров?
3. Какие операторы используются для решения систем линейных
алгебраических уравнений?
4. При каких условиях система линейных алгебраических уравнений имеет
решение?
2 3
1 2 3
1 3 1 2
а) 4 б) 4 в) 4
1
2 1
3
1 4 2 4
4
г) 3 2 3 е)
1 д)
2 4
3
ж)
Рис. 1
Таблица 1
Ветвь Граф схемы
а б в г д е ж
узлы узлы узлы узлы узлы узлы узлы
1 1 2 2 3 3 4 3 2 3 4 4 1 4 3
2 1 3 4 2 2 1 1 3 1 4 2 4 1 3
3 1 4 3 1 3 2 3 4 4 2 4 3 2 1
4 2 3 1 4 4 1 1 4 1 2 3 1 3 2
37
5 4 3 3 4 4 2 2 1 2 3 1 2 1 4
6 2 4 2 1 3 1 2 4 3 1 3 2 4 2
Таблица 2
№ R1 R2 R3 R4 R5 R6 R
V
Ом
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
5
8
3
7
2
9
10
6
7
8
6
5
8
3
5
8
4
3
2
9
6
7
8
5
8
25
23
21
32
16
18
16
24
27
17
13
19
12
21
24
23
25
29
14
15
13
18
17
15
19
20
19
12
14
16
18
17
19
15
13
11
27
21
13
18
17
18
19
22
27
21
17
25
14
12
11
9
7
9
5
7
10
11
12
16
8
14
9
11
14
7
8
9
5
9
10
9
6
7
8
19
20
21
22
23
18
15
13
12
14
16
18
19
20
6
16
14
18
12
22
21
18
16
14
15
14
13
12
11
10
11
12
13
14
15
10
9
17
15
10
11
9
10
8
19
18
21
12
21
22
100
150
120
130
90
100
110
130
120
200
160
140
180
170
100
120
200
120
190
140
150
160
170
180
190
Таблица 3
№ Е1 Е3 Е5 J1 J4 J6
B A
1
2
3
4
5
-50
25
44
28
-31
10
-12
14
-16
18
31
29
-27
25
23
1
0.5
1
-1.5
1.2
0
1
0
1
0
1
0
1
0
0.5
38
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
66
-46
54
-12
19
22
24
-28
37
34
39
-50
28
19
24
-36
38
44
27
-31
-20
19
-17
15
-13
11
-10
12
-14
16
-18
20
-21
19
-17
15
-13
11
-13
15
21
-19
20
22
24
-26
28
30
31
-19
117
21
23
-25
27
29
18
-16
34
24
-0.4
0.7
-2
1.3
-1.5
-1.1
1.3
1
-1.3
0.5
0.8
0.9
1
-1.8
1.4
-1.2
2
-2.1
1.1
-1.4
1
0
1
0
1
0
1
0
1
0
1
0
1
0
1
0
1
0
1
0
0
0.6
0
1
0
0.6
0
0.4
0
0.3
0
0.2
0
0.1
0
0.8
0
0.5
0
0.7
Литература
1. С. В. Поршнев MATLAB 7. Основы работы и программирования. Учебник
– М.: ООО <<Бином-Пресс>>, 2006г. -320 с.: ил.
2. Курбатов Е.А. MATLAB 7. Самоучитель. – М.: Издательский дом
“Вильямс”, 2006. – 256 с.: ил.
3. Половко А.М., Бутусов П.Н. MATLAB для студента. – СПб.: БХВ-
Петербург, 2005. – 320 с.: ил.
4. В.П. Дьяконов MATLAB 6.5 SP1/7 + Simulink 5/6. Основы применения.
Серия << Библиотека профессионала >>.-М.: СОЛОН-пресс, 2005. -800 с.: ил
5. Алексеев Е.Р., Чеснокова О.В. Решение задач вычислительной математики
в пакетах Mathcad 12, MATLAB 7, Maple 9 /Алексеев Е.Р., Чеснокова О.В. –
М.: НТ Пресс, 2006. – 496с. : ил.
Лабораторная работа № 1 - 5
Решение математических задач
Цель работы: изучить операции над полиномами;
дифференцированиe функции одной переменной.
Введение
39
В теории управления и теоретической электротехнике задача
определения корней и вычисление значений полиномов в некоторой точке
возникает при расчёте переходных процессов и анализе устойчивости систем
управления.
Полином или многочлен – это выражения вида
....)(
;...)(
0
1
1
1
10
axaxaxP
axaxaxP
n
n
n
n
n
nn
В MATLAB Р(х) рассматривается как вектор, элементами которого
являются коэффициенты полинома, число элементов всегда больше на
единицу степени полинома n, нулевые коэффициенты должны содержаться в
векторе. С полиномами можно выполнять операции умножения, деления,
интегрирования и дифференцирования. Рассмотрим некоторые из них на
конкретных примерах.
Пусть заданы два полинома:
.201275)(
;43)(
234
2
pppppD
pppH
Для задания полиномов в MATLAB необходимо выполнить
следующую команду:
>> H=[1 3 4]
H =
1 3 4
>> D=[1 5 7 12 20]
D =
1 5 7 12 20
С помощью оператора polyval(…) можно вычислить значение полинома
в точке р = р
i
= 5i:
>> polyval(H,5i)
ans =
-21.0000 +15.0000i
Аргумент полинома может быть задан как вектор или матрица. Для
каждого элемента вектора или матрицы будет вычислено значение полинома,
в результате получается вектор или матрица того же размера, что и аргумент
р.
>> x=[1 3 5i 4-7i]
x =
1.0000 3.0000 0 + 5.0000i 4.0000 - 7.0000i
>> polyval(D,x)
ans =
1.0e+003 *
0.0450 0.3350 0.4700 - 0.5650i -4.8300 + 3.2550i
>> s=[2 1 6;3 5 4]
s =
40