Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
Cap. 2.- Matrizes e
Sistemas Lineares
2.1. Definição
Matriz é um conjunto organizado de números dispostos em linhas e
colunas.
Representações
[
a11 a 12 ⋯ a 1n
a 2n
A= a 21
⋮
⋮
a m1 a m2 ⋯ a mn
Matriz retangular A,
m x n (eme por ene)
]
A ou [ A] ou  A ou ∥A∥
linha = rows
coluna = columns
a ij é o elemento da matriz localizado na linha i e na coluna j
A = a ij  para i = 1 m e j = 1 n
2.2. Tipos
Matriz linha
Matriz coluna
Matriz quadrada de
ordem n
A=[ 1 2 3 ]
m=1
[]
1
A= 2
3
n=1
[ ]
1 2 3
A= 4 5 6
7 8 9
m=n
Os elementos da diagonal principal são: a ij para i = j
[1 5 9 ]
Os elementos da diagonal secundária são: a ij para i + j = n + 1
[3 5 7 ]
Matriz unitária
m=n=1
Matriz diagonal
Os elementos são: a ij = 0 para i≠ j
José Eduardo Mautone Barros
A=[ 3 ]
26/10/12
[ ]
1 0 0
A= 0 5 0
0 0 9
1/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.2. Tipos(cont.)
Matriz identidade
a ij = 1 para i= j
a ij = 0 para i≠ j
É a matriz diagonal onde:
[ ]
1 0 0
I 3= 0 1 0
0 0 1
Matriz triangular
superior (U)
(“upper”)
Os elementos abaixo da diagonal principal são nulos.
1 2 3
U= 0 8 5
0 0 2
Matriz triangular
inferior (L)
(“lower”)
Os elementos acima da diagonal principal são nulos.
2 0 0
L= 3 5 0
1 2 1
Matriz nula
Todos os elementos são nulos: a ij = 0 V i e j
Matriz oposta
A = -B
Matriz idêntica
A=B
A é oposta de B se: a ij = −b ij V i e j
A é idêntica a B se: a ij = bij
Matriz cheia
Matriz esparsa
Matriz de banda
[
[
0 0
0 0
]
A=
1 −3
7 2
B=
−1
3
−7 − 2
]
a=1 ; b=2 ; c=5 ; d =7
São matrizes com a maior parte dos elementos não nulos.
São matrizes com a maior parte dos elementos nulos
São matrizes quadradas esparsas cuja diagonal principal e algumas
diagonais paralelas a principal são compostas de elementos não nulos.
[ ]
1
1
0
0
Matriz tridiagonal
José Eduardo Mautone Barros
[ ]
N=
V iej
[ ][ ]
a b
1 2
=
c d
5 7
[ ]
[ ]
26/10/12
1 0 0
2 2 0
2 3 3
0 3 4
2/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.3. Operações
Adição
As matrizes são do mesmo tamanho m x n.
c ij = aij  bij
C=A+B
V i e j
[ ][ ][ ]
1 2
4 7
5 9

=
1 3
5 8
6 11
[ 2 3 0 1 ][ 7 2 5 3 ] =[ 9 5 5 4 ]
Propriedades
A+B=B+A
A + (B + C) = (A + B) +C
A+0=A
A+(-A) = 0
Subtração
C=A-B
C = A - B = A + (-B)
Multiplicação por um
número k
C=kB
Obs: a e b podem ser
números complexos
Multiplicação
C = A.B
[ ][ ][ ]
4 7
7 1
−3 6
2 3 − 3 2 = −1 1
0 5
7 0
−7 5
c ij =k b ij
3
Propriedades
comutativa
associativa
[
][
2
3
6
9
=
1 −4
3 −12
Vi e j
[] [ ][ ]
2
1
8
2 3 4 1 = 10
5
−2
2
]
a (b A) = (a b) A
a (A + B) = a A + a B
(a +b) A = a A + b A
1.A = A
A = ( aij ) m× p
B = ( b jk ) p×n
Definição indicial
C = cik  m×n
onde cik = a i1 b1k  ai2 b2k  a i3 b3k    a ip b pk
Obs: matrizes
quadradas devem ter a
mesma ordem para
poderem ser
multiplicadas
José Eduardo Mautone Barros
p
c ik =
26/10/12
∑ aij b jk
j=1
3/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.3. Operações (cont.)
Multiplicação
C = A.B
Definição esquemática
b1k
b2k
a i1 a i2 a i3 … a ip ⋅ b3k
→
...
b pk
i−ésima linha de A k−ésima coluna de B
cik
elemento ik de C
p
c ik =
Wikipédia, 2009
Propriedades
[ ][ ] [
[ ][ ] [
∑ a ij b jk
j=1
A.B=
1 2 4 6
18 22
.
=
3 5 7 8
47 58
B.A=
4 6 1 2
22 38
.
=
7 8 3 5
31 54
A.B ≠ B.A
A.B = 0
]
]
[
][ ] [ ]
4
1 1 2 .
14
0 =
2 3 1
13
5
2 x3
3 x1 2 x 1
não é comutativa
≠ > A = 0 ou B = 0
[ ][ ] [ ]
1 0 0 0
0 0
.
=
1 0 1 1
0 0
(A.B).C = A.(B.C)
(A+B).C = A.C+B.C
C.(A+B) = C.A+C.B
(k.A).B = A.(k.B) = k.(A.B)
A.In = Im.A = A
Matriz Transposta
At
Propriedades
Matriz simétrica
Matriz anti-simétrica
associativa
distributiva a direita
distributiva a esquerda
k = constante real ou imaginária
A é uma matriz m x n
At = (bji) , tipo m x n, é a matriz transposta de A = (aij), tipo n x m onde,
[ ] [
[ ][ ] [ ]
[ ][ ] [ ]
[ ]
2 3
bji = aij V i e j A=
4 1
0 6
(A+B)t = At + Bt
(kA)t = kAt
(A.B)t = Bt.At
A.B=
At = 2 4 0
3 1 6
1 5 7
B= 4 3 2
9 6 8
] [ ]
1 4 9
Bt = 5 3 6
7 2 8
2 4 3 3
22 22
.
=
1 3 4 4
15 15
B t . At = 3 4 . 2 1 = 22 15 = A.Bt
3 4 4 3
22 15
É a matriz quadrada cuja transposta é igual a matriz original:
At = A
ou seja, aij = aji V i e j
[ ]
A=
1 2
2 3
É a matriz quadrada cuja transposta é igual a oposta da matriz original:
At = -A
ou seja, aij = -aji V i e j
A= 0 1
−1 0
[
José Eduardo Mautone Barros
26/10/12
4/30
]
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.4. Determinantes
Definição
Seja uma matriz quadrada M=(aij), de ordem n, chamamos determinante
de M, simbolizado por:
∣
∣
a 11 a12 a 13 ... a1n
a
a
a 23 ... a 2n
det M = 21 22
:
:
:
:
:
a n1 a n2 a n3 ... a nn
a um número calculado por:
a) Para n=1 então M=(a11) e det M = a11
n
det M =∑ (−1)i +1 a i1 Di1
b) Para n ≥ 2 então
i=1
onde,
Menor complementar
do elemento ai1 da
matriz M
Matriz quadrada de
ordem 2
Di1 é o determinante da matriz que se obtém de M, suprimindo-se a linha i
e a coluna 1.
∣
det M =
∣
a 11 a12
=(−1)1+1 a 11 D 11+(−1)2+1 a 21 D 21
a 21 a 22
∣ ∣
∣ ∣
Menor Complementar
D 11= : :
: a 22
Determinante 2x2
det M =a 11 a 22−a 21 a 12
D 21=
: a12
: :
Regra gráfica
∣
∣
Matriz quadrada de
ordem 3
a 11 a12 a 13
det M = a 21 a 22 a 23 =(−1)1+1 a11 D11+(−1)2 +1 a 21 D21+(−1)3 +1 a31 D31
a 31 a32 a 33
Menor Complementar
: :
:
D 11= : a 22 a 23
: a32 a 33
José Eduardo Mautone Barros
∣
∣
26/10/12
∣
∣
: a 12 a 13
D21= : :
:
: a 32 a 33
∣
∣
: a12 a 13
D 31= : a 22 a 23
: :
:
5/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.4. Determinantes
(cont.)
Determinante 3x3
det M =a 11 (a 22 a 33−a 32 a 23)−a 21 (a 12 a 33−a32 a 13)+ a 31( a 12 a 23−a 22 a 13)
det M =( a11 a 22 a 33+ a 12 a 23 a 31+ a13 a 21 a 32)−(a13 a 22 a 31+ a 11 a 23 a 32+ a12 a 21 a 33)
Regra de Sarrus
(regra gráfica)
Pierre Frédéric Sarrus
1798-1861
Matriz quadrada de
ordem n>3
Definindo Cofator do elemento aij de uma matriz quadrada de ordem n,
por
i+ j
Aij =(−1)
D ij
n
det M =∑ a i1 Ai1
Determinante nxn
Propriedades
i=1
det Mt = det M
Se uma linha ou coluna da matriz M é constituída de elementos nulos,
então det M = 0
Se multiplicarmos uma linha ou coluna da matriz M por um número k,
gerando uma nova matriz N, então det N = k det M
Se duas linhas ou colunas na matriz M forem iguais ou proporcionais,
então det M = 0
Numa matriz triangular superior ou inferior, temos
det M =a 11 a 22 a 33 ... a nn
José Eduardo Mautone Barros
26/10/12
6/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.5. Matrizes
inversíveis
Matriz Inversa (A-1)
Matriz inversa A-1 da matriz quadrada A, de ordem n, é definida por
−1
−1
A.A = A . A=I n
Obs:
Se existir inversa, então a matriz A é dita inversível ou não singular;
Se não existe inversa, então a matriz A é chamada de singular;
Se a matriz A é inversível, então a sua inversa é única.
Matriz de
Cofatores (A')
É a matriz composta pelos cofatores de cada elemento da matriz A.
A '=( Aij )
Matriz Adjunta (Ᾱ)
É a transposta da matriz de cofatores da matriz A.
̄ =( A' )t
A
Teorema
Se A é uma matriz quadrada de ordem n, tal que, D = det A ≠ 0, então a
inversa de A é,
A−1 =
1
A
D
pois,
A A=A A=D I n
Gabriel Cramer
1704-1752
José Eduardo Mautone Barros
D A−1 = A
−1
D AA = AA
DIn = AA
26/10/12
D A−1 = A
−1
D A A= AA
DIn= AA
7/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.6. Matrizes em
Planilhas e no SciLab
José Eduardo Mautone Barros
26/10/12
8/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.7. Resolução de
Sistemas Lineares
2.7.1. Representação
Notação usual
a 11 x1a 12 x 2 a1n x n=b 1
S n = a 21 x 1a 22 x 2a 2n x n=b 2
⋮ ⋮ ⋮ ⋮ ⋮ ⋮
a n1 x 1a n2 x 2a nn x n=b n
n
Notação indicial
S n =∑ a ij x j=b i
i = 1..n
j=1
Notação matricial
Matriz de coeficientes
A.X =B
[
a 11
A= a 21
⋮
a n1
Método de Cramer
José Eduardo Mautone Barros
a 1n
a 2n
⋮
a nn
]
[]
[]
x1
X = x2
⋮
xn
Matriz de incógnitas
Se D = det A ≠ 0
⋯
⋯
⋮
⋯
b1
B= b 2
⋮
bn
Matriz de termos
independentes
2.7.2. Solução
a12
a 22
⋮
a n2
A.X =B
A−1 . A.X = A−1 . B
 A−1 . A= A−1 . B
I n . X = A−1 . B
X = D1 A . B
26/10/12
9/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.7.3. Classificação
Sistema determinado det A ≠ 0
Exemplo
solução única
x2
Sistema linear 2x2
a 11 x1a 12 x 2 =b1
a 21 x 1a 22 x 2=b 2
x1
Sistema indeterminado det A = 0
múltiplas soluções
x2
Uma equação
inconsistente!
x1
Sistema incompatível
sem solução
x2
Uma equação
redundante!
x1
Sistema homogêneo
B=(0)
Solução trivial
x2
Um sistema linear
homogêneo sempre
admite a solução
trivial!
José Eduardo Mautone Barros
X=(0)
x1
26/10/12
10/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.7.3. Classificação
(cont.)
Sistema triangular
Admite um método de solução simplificado.
O sistema triangular superior é resolvido por substituição retroativa.
Exemplo
A.X =B
[
3
0
0
0
][ ] [ ]
x1
4 −5 1
−10
1 1 −2 . x 2 = −1
0 4 −5 x 3
3
0 0
2
2
x4
x 4 =2 /2=1
4x 3−5⋅1=3⇒ x 3=2
x 2 2−2⋅1=−1⇒ x 2=−1
3x 14⋅−1−5⋅21=−10 ⇒ x1 =1
X t =[ 1 −1
2
1]
Obs
O sistema triangular inferior é resolvido por substituição progressiva.
2.7.4. Operações com
Sistemas Lineares
Realizando as seguintes operações sobre as linhas de um sistema linear
obtém-se uma sistema linear equivalente ao primeiro, ou seja, um novo
sistema que possui a mesma solução:
•
•
•
Obs
2.7.5. Aplicações
A troca de colunas na matriz de coeficientes somente altera a ordem dos
termos na matriz solução X.
•
•
•
•
•
•
•
•
José Eduardo Mautone Barros
Trocar a posição de duas equações do sistema;
Multiplicar uma equação do sistema por uma constante não nula;
Adicionar duas equações do sistema.
Elétrica e Eletrônica
Economia
Engenharia Civil
Engenharia Térmica
Engenharia Aeronáutica
Engenharia de Estruturas
Otimização
Computação Gráfica
26/10/12
11/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.7.5. Aplicações
(cont.)
Modelo de um sistema mecânico
F1
P1
F1
x=0
K1
K1
x1
m1
F2
a
2 Lei de Newton no equilíbrio
K2
K2
P1
F2
W1
P2
F2
F2
x2
m2
∑ Fext =0
K3
K3
P2
x3
F3
F3
P3
F3
F3
m3
P3
Diagrama de corpo livre
W2
Massa 1:
∑ F x=2F 1−P 1−2F 2 −W 1
Massa 2:
∑ F x =2F 2−P 2−2F 3−W 2
Massa 3:
∑ F x=2F 3−P 3−W 3
W3
F 1=K 1 x 1
F 2=K 2  x 2−x 1 
F 3=K 3  x 3−x 2
Força de uma mola
K é a constante da
mola
Força devido ao
Peso Próprio
g é aceleração da
gravidade
W 1=m1 g
W 2=m2 g
W 3=m3 g
O sistema de equações lineares que modelam o problema é:
2  K 1 K 2 x 1−2K 2 x 2= P1 m1 g
−2K 2 x1 2  K 2K 3 x 2−2K 3 x 3= P 2m2 g
−2K 3 x 2 2K 3 X 3=P 3m3 g
A matriz de
coeficientes K é
chamada de matriz de
rigidez do sistema e é
simétrica.
2
[
][ ] [ ] [ ]
K 1K 2 −K 2
0
x1
P1
m1 g
=

−K 2
K 2 K 3 −K 3 x 2
P2
m2 g
0
−K 3
K 3 x3
P3
m3 g
2 K X =P W
O sistema homogêneo
P=0 é usado em
análise de vibrações
naturais.
José Eduardo Mautone Barros
●
●
●
para P = 0 o sistema está em equilíbrio devido ao peso próprio;
para uma dada carga P os deslocamentos X são únicos;
para uma carga cíclica P os deslocamentos X também são cíclicos.
26/10/12
12/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.7.5. Aplicações
(cont.)
Ponte de Wheatstone
Instrumentação:
Uma Ponte de Wheatstone (ao lado) é um circuito elétrico usado para
medição de sinais de vários tipos de sensores (de termistores a células de
carga). O medidor de tensão, de resistência Rg, fica entre os terminais D e
B da ponte. A fonte fornece uma força eletromotriz E ao circuito.
Vamos mostrar que para uma condição de equilíbrio na ponte, temos:
I g =0
⇒
R4 R 1
=
R3 R 2
As leis que regem o fenômeno físico são:
Ao longo de qualquer
circuito envolvendo a
fonte.
Lei de Ohm
V =R.I
Nós do circuito =
pontos A, B, C e D na
figura
Lei de Kirchhoff
∑ I nó=0
O sistema de equações lineares que modelam o sistema é:
I 1−I 2− I g =0
I 3−I 4I g =0
R1 I 1R 2 I 2= E
R 3 I 3 R4 I 4=E
R1 I 1R3 I 3R g I g =E
Nó D
Nó B
Circuito ADC
Circuito ABC
Circuito ADBC
Na notação matricial:
[
][ ] [ ]
1 −1 0
0 −1
0
0 1 −1 1
R1 R2 0
0
0
0
0 R 3 R4 0
R1 0 R 3 0 R g
ou seja,
Usar substituição
progressiva!
I1
0
I2
0
I3 = E
E
I4
E
Ig
A X =B
Resolvendo o sistema para Ig = 0 , obtemos,
R4 R1
=
R3 R2
José Eduardo Mautone Barros
26/10/12
13/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.7.6. Métodos
Numéricos de Solução
de Sistemas Lineares
Métodos Diretos
A solução é encontrada por métodos algébricos com um número fixos de
operações. A solução é exata.
Recomendados para:
•
•
Exemplos de métodos
diretos
•
•
•
•
Métodos Indiretos
Sistemas lineares pequenos (n<=1000).
Matriz de coeficientes do tipo matriz cheia, onde a maioria dos
elementos são não nulos (aij ≠ 0).
Método da Eliminação de Gauss
Método da Eliminação de Gauss-Jordan
Método da Inversão da Matriz de Coeficientes
Método da Decomposição LU
A solução é encontrada por tentativa e erro, através de um processo
iterativo. Uma solução é assumida e substituída no sistema de equações
para o cálculo do erro. Este erro é usado para melhorar a estimativa. O
procedimento é repetido até que o erro calculado seja menor que um valor
pré-definido. A solução final é aproximada.
Recomendados para:
•
•
Exemplos de métodos
indiretos
•
•
•
(SOR = successive
over relaxation)
2.7.6.1.. Métodos
Diretos
José Eduardo Mautone Barros
•
Sistemas lineares grandes (n>1000).
Matriz de coeficientes do tipo matriz esparsa , onde a maioria dos
elementos são nulos (aij =0).
Método de Iteração de Jacobi
Método de Iteração de Gauss-Siedel
Método da Relaxação
Método da Super-Relaxação Sucessiva (SOR)
Serão discutidos os Métodos de Gauss e o da Decomposição LU.
26/10/12
14/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.7.6.11. Método de
Gauss
Método da Eliminação
de Gauss
É um método direto de solução de sistemas lineares.
Consiste em transformar um sistema de equações lineares em um outro
sistema triangular superior, equivalente ao primeiro, e de solução direta
por substituição retroativa.
Método da
Triangularização de
Gauss
[
a 11
a 21
⋮
a n1
a12
a 22
⋮
⋯
⋯
⋯
⋮
⋯
a1n
a 2n
⋮
a nn
b1
b2
⋮
bn
] [
→
1
0
⋮
0
c12
1
⋮
0
⋯
⋯
⋮
⋯
Etapas:
1) Escrever a matriz aumentada
Carl Friedrich Gauss
1777-1855
Fórmula generalizada
de transformação:
•
•
José Eduardo Mautone Barros
⋯
⋯
⋮
⋯
c 1n
c 2n
⋮
c nn
c1
c2
⋮
cn
m021=−c 021 /c 011
m031=−c 031 /c011
3) Transformar as linhas L para obter a nova matriz aumentada C1
L11 =L01
L12=m 021 L01 L02
L13=m 031 L 01 L03
4) Repetir as operações de pivotamento e transformação para o novo pivô
até chegar a última linha da matriz aumentada.
m132=−c 132 /c 122
Obs.: nenhum
elemento da diagonal
principal pode ser
nulo.
c12
c22
⋮
⋯
]
Escolher o pivô (elemento da diagonal principal)
c011
Encontrar os multiplicadores para eliminar os termos c021 e c031
Lik =mijk−1 Lk−1
+ L k−1
j
i
para i = k+1 .. n
para j = k
k = índice de iteração;
i = índice da linha;
j = índice da coluna.
[
c 11
c 21
⋮
c n1
d1
d2
⋮
dn
2) Pivotamento
k−1
ij
k−1
jj
−c
mkij −1=
c
para i = k+1 .. n
para j = k
C0= [A : B]
c 1n
c 2n
⋮
1
L21= L11
L22= L12
L23=m132 L 12+ L 13
5) Resolver o sistema triangular superior por substituição retroativa.
26/10/12
15/30
]
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.7.6.1.1. Método de
Gauss (cont.)
Exemplo
Sistema de ordem 3
AX=B
2x 1 3x2 −x 3=5
4x 14x 2 −3x 3=3
2x 1 −3x2 x 3=−1
[
[
]
[
]
C = [A : B]
2
3 −1 ⋮ 5
C 0= 4
4 −3 ⋮ 3
2 −3 1 ⋮ −1
L12 =
-2*[ 2 3 -1 : 5] +
[ 4 4 -3 : 3 ]
= [ 0 -2 -1 : -7]
2
3 −1 ⋮ 5
C 1 = 0 −2 −1 ⋮ −7
0 −6
2 ⋮ −6
Matriz triangular
superior
2
3 −1 ⋮ 5
C 2 = 0 −2 −1 ⋮ −7
0
0
5 ⋮ 15
Solução por
substituição retroativa
5x 3=15
−2x 2− x 3=−7
2x 13x 2−x 3=5
0
Solução Final
[]
1
X= 2
3
Obs.:
José Eduardo Mautone Barros
ou
]
m021=−c021 /c011=−2
m 031=−c031 /c011=−1
L11 =L01
L 12 =m021 L01 + L 02
L13 =m031 L 01 + L 03
L21 =L 11
L 22 =L 12
m132=−c132 /c122=−3 L 23 =m 132 L12 + L 13
X t =[ 1 2 3 ]
n = 3 exige 31 operações aritméticas
26/10/12
16/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.7.6.1.1. Método de
Gauss (cont.)
Exemplo
Matriz triangular
superior
Solução por
substituição retroativa
Sistema de ordem 4
[
[
[
[
[]
3
2
0
1
0
9
8 −3
4
C =
−6
4 −8
0
3 −8
3 −4
⋮
3
⋮
6
⋮ −16
⋮ 18
3
2
0
1
1
0
2 −3
1
C =
0
8 −8
2
0 −10
3 −5
⋮ 3
⋮ −3
⋮ 10
⋮ 15
3
2
0
C =
0
0
2
0
1 ⋮ 3
2 −3
1 ⋮ −3
0
4 −2 ⋮ 2
0 −12
0 ⋮ 0
3
0
C =
0
0
2 0
1 ⋮ 3
2 −3
1 ⋮ −3
0
4 −2 ⋮ 2
0
0 −6 ⋮ 6
3
2
−1
X=
0
−1
Obs.:
ou
t
X =[ 2 −1
]
]
]
m 021=−c021/ c 011=−3
m 031=−c031 /c011= 2
0
0
0
m41 =−c 41/ c 11=−1
L11 =L 01
L12 =m 021 L01 + L02
L13 =m 031 L01 + L 03
1
0
0
0
L4 =m41 L1 + L 4
L21 =L 11
L 22 =L12
m 132=−c 132 /c 122=−4
1
1
1
m42 =−c 42 /c 22= 5
L 23 =m132 L12 + L 13
2
1
1
1
L 4 =m42 L2 + L 4
L31 =L 21
L32 =L 22
L33 =L 23
0
2
2
3
2
2
2
m43=−c 43 /c 33= 3 L 4 =m43 L3 + L 4
]
0 −1 ]
n = 4 exige 76 operações aritméticas
n = 5 exige 145 operações aritméticas
Exercício: Obter uma equação para calcular o número de operações
aritméticas necessárias para a solução de um sistema linear de ordem n
pelo método de eliminação de Gauss.
José Eduardo Mautone Barros
26/10/12
17/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.7.6.1.1. Método de
Gauss (cont.)
Algorítimo de Triangularização de Gauss para Sistema Lineares
No SciLab:
function x = GaussElim(n,a,b)
// Matriz aumentada
Entradas:
Ordem do sistema linear
Matriz de coeficientes
Matriz de termos independentes
n
a[n,n]
b[n,1]
Saída:
Matriz de incógnitas
x[n,1]
Início
c = [a b];
// Definir os termos da matriz aumentada
c[n,n+1]= a[n,n]:b[n,1];
// Triangularização da matriz
// aumentada
// Triangularização da Matriz aumentada
for k=1:n-1
for i=k+1:n
mik=c(i,k)/c(k,k);
c(i,k) = 0;
for j = k+1:n+1
c(i,j)=c(i,j)-mik*c(k,j);
end
end
end
// Substituição retroativa
x=zeros(n,1);
x(n)=c(n,n+1)/c(n,n);
for i=n-1:-1:1
soma = 0;
for j=i+1:n
soma = soma +c(i,j)*x(j);
end
x(i)=(c(i,n+1)-soma)/c(i,i);
end
Para k=1 até n-1 faça
início
Para i=k+1 até n faça
início
mik=c(i,k)/c(k,k);
c(i,k) = 0;
Para j = k+1 até n+1 faça
início
c(i,j)=c(i,j)-mik*c(k,j);
fim;
fim;
fim;
// Substituição retroativa
endfunction
x(n)=c(n,n+1)/c(n,n);
Para i=n-1 até 1 faça
início
soma = 0;
Para j=i+1 até n faça
início
soma = soma +c(i,j)*x(j);
fim;
x(i)=(c(i,n+1)-soma)/c(i,i);
fim;
Mostre a matriz x[n,1];
fim.
José Eduardo Mautone Barros
26/10/12
18/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.7.6.1.2. Método de
Decomposição LU
AX =B
[
2
A= 1
0
3 −1
0 2
3 −1
] [ ] []
x1
X = x2
x3
4
B= 3
2
A=LU
LUX =B
LY =B
UX =Y
A matriz U é única!
LY =B Y  UX =Y  X
Método de solução
l ij =1 se i= j
l ij =0 se i j
u ij =0 se i j
Decomposição usando
a igualdade LU = A
A= LU
[
1 0 0
L= l 21 1 0
l 31 l 32 1
] [
u 11 u 12 u 13
U= 0
u 22 u 23
0
0
u 33
] [
a 11 a 12 a 13
A= a 21 a 22 a 23
a 31 a 32 a 33
u 11=a 11 u 12=a12 u 13 =a 13
l 21 u 11=a 21 ⇒ l 21 =a 21 / u11
l 21 u 12 +u 22=a 22 ⇒ u 22=a 22−l 21 u 12
l 21 u 13+ u 23=a 23 ⇒ u 23=a 23−l 21 u 13
l 31 u 11=a 31 ⇒ l 31 =a 31 /u11
l 31 u 12+ l 32 u 22=a 32 ⇒ l 32 =(a32−l 31 u12 )/u 22
l 31 u 13+ l 32 u 23+ u33=a 33 ⇒ u 33=a 33−l 31 u 13−l 32 u 23
Forma indicial
Obs: na fórmula de uij
j-1<i significa
ignorar os termos do
somatório onde k≥i
]
j−1< i
para i⩽ j
para i> j
u ij =a ij se j=1
e
a ij
se j=1
u jj
e
l ij =
u ij =a ij −
∑
k=1
j−1
l ik u kj se j> 1
l ij =(aij −∑ l ik u kj) /u jj se j> 1
k=1
LY =B
Substituição
progressiva
Substituição
retroativa
[
1
0 0
L= 1/2 1 0
0
−2 1
4
B= 3
2
UX =Y
[
2 3
−1
U = 0 −3/2 5/2
0 0
4
] [ ] []
x1
X = x2
x3
4
Y= 1
4
X t =[ 1 1 1 ]
Solução
José Eduardo Mautone Barros
] [ ] []
y1
Y = y2
y3
26/10/12
19/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.7.6.1.2. Método de
Decomposição LU
(cont.)
Rotina de
Decomposição LU em
código do SciLab
//Rotina de Decomposição LU para Sistema Lineares
//entrada Ordem do sistema linear
n
//
Matriz de coeficientes
a[n,n]
//
Matriz de termos independentes b[n,1]
//saída
Matriz de incógnitas
x[n,1]
function x = DecoLU(n,a,b)
// Decomposição de A em L(matriz triangular inferior) e U(matriz triangular superior)
// LU=A
l = zeros(n,n);
u = zeros(n,n);
for i=1:n
l(i,i)=1;
end
// zerar matrizes L e U
// diagonal de L igual a 1
j=1;
// cálculo dos elementos de L e U para j=1
for i=1:n
if i<=j then
u(i,j)=a(i,j);
else
l(i,j)=a(i,j)/u(j,j);
end
end
Obs: O programa para
o cálculo dos uij ignora
os termos de k≥i
automaticamente, pois
estes ainda são nulos
até o seu cálculo.
for i=1:n
// cálculo dos elementos de L e U para j>1
for j=2:n
SumLU=0;
for k=1:j-1
SumLU=SumLU+l(i,k)*u(k,j);
end
if i<=j then
u(i,j)=a(i,j)-SumLU;
else
l(i,j)=(a(i,j)-SumLU)/u(j,j);
end
end
end
// Substituição progressiva LY=B
y=zeros(n,1);
y(1)=b(1);
for i=2:n
SumLY=0;
for j=1:i-1
SumLY = SumLY + l(i,j)*y(j);
end
y(i)=b(i)-SumLY;
end
// Substituição retroativa UX=Y
x=zeros(n,1);
x(n)=y(n)/u(n,n);
for i=n-1:-1:1
SumUX = 0;
for j=i+1:n
SumUX = SumUX +u(i,j)*x(j);
end
x(i)=(y(i)-SumUX)/u(i,i);
end
endfunction
// Jose Eduardo Mautone Barros 22/04/2010
José Eduardo Mautone Barros
26/10/12
20/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.7.6.2. Métodos
Indiretos (Iterativos)
Exemplo:
Os métodos iterativos consistem em transformar o sistema de equações
lineares original para uma outra forma que permita obter novas
estimativas de valores do vetor de incógnitas X a partir de uma estimativa
anterior de valores do vetor X.
AX =B
AX −B=0
AX IX −B= IX
X = AI  X −B
A X =B
para
X =F X D
A partir de uma aproximação inicial:
X  0 t =[ x 01  x20 x30  ... x 0n ]
obtemos a nova estimativa ,
X 1=F X 0 D
e repete-se até que,
máx∣x ki 1− x ki ∣≤
ou
kM
onde,
2.7.6.2.1. Método de
Jacobi
ε = tolerância na solução
M = número máximo de iterações
Seja o sistema de equações lineares (LES),
a 11 x1a 12 x 2 a 13 x 3...a1n x n=b 1
a 21 x 1a 22 x 2a 23 x 3...a 2n x n=b2
...
a n1 x 1a n2 x 2a n3 x 3...a nn x n=bn
explicita-se as incógnitas x da seguinte forma:
Carl Gustav Jakob
Jacobi 1804-1851
Obs: aii ≠ 0 V i
Senão é necessário
reagrupar as equações
do sistema original.
José Eduardo Mautone Barros
b 1−a 12 x 2 a 13 x 3...a 1n x n 
a 11
b − a 21 x 1a 23 x 3...a 2n x n 
x2 = 2
a 22
...
b − a n1 x 1an2 x 2 ...a n n−1 x n −1 
xn = n
a nn
x 1=
26/10/12
21/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.7.6.2.1. Método de
Jacobi (cont.)
O método iterativo de Jacobi consiste em:
a) Partindo-se de uma aproximação inicial
X 0= x 01 , x 02 , x03 , , x 0n t
b) Calcula-se a sequência de aproximações
X1,X2, X3, ..., Xk
utilizando as equações:
1
b −a x k −a x k −a x k −−a 1n x kn 
a 11 1 12 2 13 3 14 4
1
x2k1 = b 2−a 21 x k1 −a 23 x k3−a 2 x k4−−a 2n x kn 
a 22
1
x3k1 = b3−a 31 x k1 −a 32 x k2−a 34 x k4−−a 3n x kn 
a 33
          
1
xnk1 = b n−a n1 x k1 −a n2 x k2−−a n n −1  x nk −1 
a nn
x1k1 =
c) Continuar a gerar aproximações até que uma das seguintes condições
for satisfeita:
k 1
k
máx∣x i − x i ∣≤
ou
kM
onde,
ε = tolerância
M = número máximo de iterações
Método do Resíduo
Ri(k)
 k1
xi
É mais atual!
k
=x i 
Rki 
a ii
i=1..n
n
Ri k=bi−∑ a ij x jk
i=1..n
j =1
José Eduardo Mautone Barros
26/10/12
22/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.7.6.2.2. Método de
Gauss-Siedel
Seja o sistema:
AX=B
O método iterativo de Gauss-Siedel consiste em:
a) Partindo-se de uma aproximação inicial
0
0
0
0
0 t
X = x 1 , x 2 , x 3 , , x n 
b) Calcula-se a sequência de aproximações
X1,X2, X3, ..., Xk
utilizando as equações:
1
k
k
k
k
b −a x −a x −a x −−a 1n x n 
a 11 1 12 2 13 3 14 4
1
x2k1 = b 2−a 21 x k1 1 −a 23 x3k −a 2 x4k −−a 2n x kn 
a 22
1
k
k
x3k1 = b3−a 31 x k1 1−a 32 x k1
2 −a34 x 4 −−a 3n x n 
a 33
          
1
k1
k 1
xnk1 = b n−a n1 x k1
1 −a n2 x 2 −−a n n−1 x n−1 
a nn
 k1
x1
=
c) Continuar a gerar aproximações até que uma das seguintes condições
for satisfeita:
máx∣x ki 1− x ki ∣≤
ou
kM
onde,
ε = tolerância
M = número máximo de iterações
Método do Resíduo
Ri(k)
É mais atual!
 k1
xi
k
=x i 
Rki 
a ii
i −1
n
j =1
j =i
Ri k=bi−∑ a ij x jk1−∑ a ij x jk
José Eduardo Mautone Barros
26/10/12
i=1..n
i=1..n
23/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.7.6.2.3. Superrelaxação Sucessiva
(SOR)
Alteração do método de Gauss-Siedel para acelerar a convergência.
SOR = successive
over-relaxation
 k
x
Método do Resíduo
Ri(k)
R
= x  i
a ii
 k1
i
k 
i
i −1
R =bi−∑ a ij x
 k
i
j =1
 k1
j
i=1..n
n
−∑ a ij x jk
i=1..n
j =i
Fator de relaxação (ω)
w< 1 sub-relaxado
w> 1 super-relaxado
0<w <2
Exemplo
Método de Jacobi
2 x 1− x 2=1
x 1 +2 x 2=3
x k1 +1=(1+ x k2 )/2
x k2 +1=(3− x k1 )/2
2 x 1− x 2=1
x 1 +2 x 2=3
K
0
1
2
3
4
5
6
7
8
9
X1
0,000
0,500
1,250
1,125
0,938
0,969
1,016
1,008
0,996
0,998
E1
X2
0,000
1,500
1,250
0,875
0,938
1,031
1,016
0,992
0,996
1,002
0,500
0,750
0,125
0,188
0,031
0,047
0,008
0,012
0,002
E2
Teste
1,500
0,250
0,375
0,063
0,094
0,016
0,023
0,004
0,006
não
não
não
não
não
não
não
não
convergiu
Convergência de x1 - Jacobi
0,800
0,700
0,600
0,500
0,400
0,300
0,200
0,100
0,000
1
José Eduardo Mautone Barros
2
3
26/10/12
4
5
6
7
8
9
24/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.7.6.2. Métodos
Indiretos (cont.)
Método de
Gauss-Siedel
2 x 1− x 2=1
x 1 +2 x 2=3
x k1 +1=(1+ x k2 )/2
1
x k2 +1=(3− x k+
1 )/2
K
0
1
2
3
4
5
X1
0,000
0,500
1,125
0,969
1,008
0,998
E1
X2
0,000
1,250
0,938
1,016
0,996
1,001
0,500
0,625
0,156
0,039
0,010
E2
Teste
1,250
0,313
0,078
0,020
0,005
não
não
não
não
convergiu
Convergência de x1 - Gauss-Siedel
0,700
0,600
0,500
0,400
0,300
0,200
0,100
0,000
1
Método SOR
2 x 1− x 2=1
x 1 +2 x 2=3
Rk1 =1−2x k1 + x k2
Rk2 =3−x k1 +1−2 x k2
x k1 +1=x 1k +ω R1k /2
x k2 +1=x k2 +ω R2k /2
2
Ômega =
K
0
1
2
3
4
3
0,91
X1
0,000
0,455
1,023
1,004
1,000
4
R1
X2
0,000
1,158
1,004
0,999
1,000
1,000
1,248
-0,042
-0,009
5
R2
Teste
2,545
-0,339
-0,011
0,003
não
não
não
convergiu
Convergência de x1 - SOR
1,400
1,200
1,000
0,800
0,600
0,400
0,200
0,000
-0,200
1
José Eduardo Mautone Barros
2
26/10/12
3
4
5
25/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.7.6.2.4. Métodos
Indiretos Implementação em
Planilhas
Método de Jacobi
•
Exemplo
 4 x1 − 2 x2 + x3 = 3

 x1 − 4 x2 + x3 = − 2
 x + 2x + 4x = 7
2
3
 1
•
Erro admitido
de 1x10-2
•
•
•
•
Método de
Gauss-Siedel
 4 x1 − 2 x2 + x3 = 3

 x1 − 4 x2 + x3 = − 2
 x + 2x + 4x = 7
2
3
 1
•
Método SOR
Omega de 1,1
(colocado numa célula
fixa $I$1)
•
•
B
X1
0
FM1
C
E1
D
X2
0
FM3
FM2
E
E2
F
X3
0
FM5
FM4
G
E3
H
Teste
FM6
FM7
Definir as seguintes fórmulas:
FM1 = 1/4*(3 + 2*D2 - F2)
FM2 = ABS(B3 - B2)
FM3 = 1/4*(2 + B2 + F2)
FM4 = ABS(D3 - D2)
FM5 = 1/4*(7 - B2 - 2*D2)
FM6 = ABS(F3 – F2)
FM7 = SE((C3>0,01) E (E3>0,01) E
(G3>0,01);”convergiu”;”não”)
Selecionar as células A2 a A3;
Estender as definições destas células até o número de iterações
desejadas;
Selecionar as células B3 a H3;
Estender as definições destas células até o número de iterações
desejadas.
Definir as seguintes fórmulas:
FM1 = 1/4*(3 + 2*D2 - F2)
FM2 = ABS(B3 - B2)
FM3 = 1/4*(2 + B3 + F2)
FM4 = ABS(D3 - D2)
FM5 = 1/4*(7 - B3 - 2*D3)
FM6 = ABS(F3 – F2)
FM7 = SE((C3>0,01) E (E3>0,01) E
(G3>0,01);”convergiu”;”não”)
Selecionar e estender as coluna A e a linha 3, colunas B a H, do
mesmo modo anterior.
A
k
0
1
2
1
2
3
4
5
 4 x1 − 2 x2 + x3 = 3

 x1 − 4 x2 + x3 = − 2
 x + 2x + 4x = 7
2
3
 1
José Eduardo Mautone Barros
A
k
0
1
2
1
2
3
4
5
B
X1
0
FM1
C
R1
FM2
D
X2
0
FM3
E
R2
FM4
F
X3
0
FM5
G
R3
H
Teste
FM6
FM7
Definir as seguintes fórmulas:
FM1 = B2+$I$1*C2/4
FM2 = 3 - 4*B2+ 2*D2 - F2
FM3 = D2+$I$1*E2/(-4)
FM4 = -2 - B3 + 4*D2 - F2
FM5 = F2+$I$1*G2/4
FM6 = 7 - B3 – 2*D3 - 4*F2
FM7 = SE((C3>0,01) E (E3>0,01) E
(G3>0,01);”convergiu”;”não”)
Selecionar e estender as coluna A e a linha 3, colunas B a H, do
mesmo modo anterior.
26/10/12
26/30
I
1,1
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.8. Problemas de
Autovalor
Seja o sistema de equações lineares:
x2
AX=B
●
x2
●
x1
Se det A ≠ 0
x1
○ Ele admite solução única
Se: det A = 0
○ Ele pode não admitir solução
○ Ele pode admitir um número infinito de soluções
○ Ele pode admitir ao menos a solução trivial, X = 0, se o
sistema for homogêneo, A X= 0
Para sistemas homogêneos, com det A = 0, só existe a solução trivial se
os coeficientes aij forem fixos. Se alguns destes coeficientes for função de
uma variável, tal como lambda (λ), existem outras soluções diferentes da
trivial. Neste caso, a matriz X é chamada de autovetor e os valores de
lambda são chamados de autovalores do sistema de equações lineares.
K
K
Exemplo
M
M
x1
x2
x
O sistema da figura acima é constituído de duas massas e duas molas
idênticas. As equações diferenciais que descrevem o seu estado (posições
das massas) são obtidas a partir do balanço de forças em cada massa do
sistema. Do seguinte modo,
2
d x1
= K  x 2−x 1 −K x 1=K x 2−2 K x 1
2
dt
d 2 x2
M
=−K  x 2−x 1 =K x 1− K x 2
2
dt
M
As solução isoladas dos sistemas massa-mola são:
ϖ = freqüência
natural
ϕ = ângulo de fase
Xi = amplitude da
oscilação da massa i
José Eduardo Mautone Barros
x 1= X 1 sen   t  
x˙1= X 1  cos   t  
2
x¨1=− X 1  sen  t  
26/10/12
x 2 = X 2 sen   t 
x˙2 = X 2  cos  t  
x¨2 =−X 2  2 sen   t 
27/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.8. Problemas de
Autovalor (cont.)
Substituindo as soluções nas equações diferenciais, temos,
 2−  X 1− X 2=0
−X 11−  X 2=0
onde, lambda foi definido como,
=
Definição arbitrária
2M
K
Na forma matricial,
[
][ ]
X1
2− 
−1
=0
−1
1−  X 2
ou,
Problema de
Autovalor
 A− I  X =0
Esta é a forma clássica do Problema de Autovalor.
A equação característica do problema de autovalor é obtida pela
condição de múltiplas soluções, não triviais, para o sistema de equações
lineares.
Equação
característica
det  A− I =0
A solução gera um polinômio da mesma ordem do número de equações
do sistema e cujas raízes são os autovalores do sistema de equações
lineares homogêneo.
Solução Numérica do
Exemplo
∣
∣
2−  −1 = 0
−1
1− 
 2−  1− −1=0
2−3  1=0
Assim,
Autovalores
José Eduardo Mautone Barros
1=2,6180
26/10/12
2=0,3820
28/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
2.8. Problemas de
Autovalor (cont.)
Para calcular os autovetores, usamos cada um dos autovalores,
Para λ1 = 2,6180
[
−0,6180
−1
−1
−1,6180
][ ]
X1
=0
X2
X 2=−0,6180 X 1
Para λ2 = 0,3820
[
1,6180
−1
−1
0,6180
][ ]
X1
=0
X2
X 2=1,6180 X 1
Considerando X1 = 1 , temos,
Autovetores
Significado físico
X1 = 1 e X2 = -0,6180 para λ1
X1= 1 e X2 = 1,6180 para λ2
●
Para λ1 as massas estão se movendo em direções opostas e a
amplitude da oscilação da segunda massa é 61,8 % da amplitude
de oscilação da primeira massa. A freqüência de oscilação é dada
por:
 1=
●
Para λ2 as massas estão se movendo na mesma direção e a
amplitude da oscilação da segunda massa é 161,8 % da amplitude
de oscilação da primeira massa. A freqüência de oscilação é dada
por:
 2=
José Eduardo Mautone Barros

1 K
M
26/10/12

2 K
M
29/30
Métodos Numéricos Aplicados à Engenharia Mecânica - EMA-084N
Métodos de solução
O problema é achar os autovalores quando o polinômio gerado é de
elevado grau, por exemplo, de ordem 300. Ou seja, possui 300 raízes.
A solução do problema de achar os autovalores pode ser feita através dos
seguintes métodos:
●
Diretos – solução pela definição ou usando uma modificação da
equação característica;
●
Indiretos – solução iterativa ou outro método de busca de raízes,
tais como,
○ Método da potência
○ Método do inverso da potência
○ Método do deslocamento de autovalores.
Solução do problema
de autovalor em
programas simbólicos
Matlab
eig(A,X)
Scilab
[erots,X] = spec(A)
Obs
Existem problemas de autovalor que não são lineares:
[ A−B   ] X =0
det [ A−B  ]=0
José Eduardo Mautone Barros
26/10/12
30/30