RESOLUÇÕES DE SISTEMAS DE EQUAÇÕES LINEARES

A resolução de sistemas lineares é um problema que surge nas mais diversas áreas. Os métodos de solução de sistemas de equações lineares podem ser divididos em dois grupos . O primeiro é o grupo dos métodos exatos ou diretos, nos quais o resultado se obtém após um número finito e determinado de operações aritméticas. O segundo grupo reúne os métodos aproximados de solução fundamentados em algoritmos iterativos.

1-SISTEMAS DE EQUAÇÕES LINEARES

Uma equação é linear se cada termo contém somente uma variável e cada variável aparece na primeira potência. A um conjunto de equações lineares se dá o nome de sistema equações lineares.

Um sistema linear de n equações e n variáveis é escrito usualmente na forma:

a11x1 + a12x2 + ... + a1n xn = b1

a21x1 + a22x2 + ...+ a2n xn = b2

........................................................................

an1x1 + an2x2 + ...+ ann xn = bn

E é representado na forma matricial por :

a11 a12 + ... + a1n x1 b1

a21 a22 + ... + a2n x2 b2

... ... .... ... ... = .....

an1 an2 + ... + ann xn bn

ou simplesmente por Ax =b, onde A é chamada de matriz dos coeficientes, b é o vetor do termo independente e x é o vetor solução.Se o conjunto ordenado de números ( α1, α2, α3 , ..., αn) satisfazer todas as equações do sistema, será denominado solução do sistema linear.Mesmo no caso mais geralem que o sitema linear envolve n equações e n variáveis, apenas uma entre as situações abaixo poderá ocorrer:

a) O sistema linear tem solução única, (sistema possível e determinado);b) O sistema linear admite infinitas soluções, (sistema possível e indeterminado) ;c) O sistema linear não admite solução, (sistema impossível ).

2- MÉTODOS DIRETOS

2.1 –MÉTODO DE GAUSS- JORDAN: o método de Gauss-Jordan exige que a matriz dos coeficientes das variáveis devem ser transformada em matriz identidade.

2x1 + x2 + 3x3 = 8

Exemplo 1: Dado o sistema 4x1 + 2x2 + 2x3 = 4 , resolver pelo método de Gauss-Jordan.

2x1 + 5x2 + 3x3 = -12

Solução

2 1 3 8 ½ L1 1 1/2 3/2 4 1 ½ 3/2 4

4 2 2 4 4 2 2 4 L2 -4L1 0 0 -4 -12

2 5 3 -12 2 5 3 -12 L3-2L1 0 4 0 -20

permutando L2 por L3 temos:

1 1/2 3/2 4 1 ½ 3/2 4 L1 -1/2L2

0 -4 0 -20 ¼ L2 0 1 0 -5

0 0 -4 -12 0 0 -4 -12 -1/4 L3

1 0 3/2 13/2 L1-3/2L3 1 0 0 2

0 1 0 -5 0 1 0 .5

0 0 1 3 0 0 1 3

Assim, o sistema de equações lineares se transforma no sistema equivalente:

x1 + 0x2 + 0x3 = 2

0x1 + x2 + 0x3 = -5 x1 =2, x2 =-5 e x3 = 3 Sistema é possível e determinado.

0x1 + 0x2 + x3 = 3

3x1 + 2 x2 - 5x3 = 8

Exemplo 2: Dado o sistema 2x1 - 4x2 - 2x3 = - 4 , resolver pelo método de Gauss-Jordan.

x1 - 2x2 - 3x3 = -4

2.2- FATORAÇÃO LU: Seja o sistema linear Ax = b.O processo de fatoração para resolução deste sistema consiste em decompor a matriz A dos coeficientes em um produto de dois ou mais fatores e, em seguida, resolver uma seqüência s]de sistemas lineares que nos conduzirá “a solução do sistema linear original. Por exemplo, se pudermos realizar a fatoração: A = CD, o sistema linear Ax = b pode ser escrito:

(CD)x = b

Se y= Dx, então resolver o sistema linear Ax = b é equivalente a resolver o sistema linear Cy = be, em seguida, o sistema linear Dx = y.

A vantagem dos processosde fatoração é que podemos resolver qualquer sistema linear que tenha A como matriz dos coeficientes. Se o vetor b for alterado, a resolução do novo sistema linear será quase que imediata.

A fatoração LU é um dos processos de fatoração mais empregados. Nesta fatoração a matriz L é triangular inferior com diagonal unitária e a matriz U é triangular superior.

Exemplo 3: Resolver o sistema linear a seguir usando a fatoração LU:

3x1 + 2 x2 + 4x3 = 1 3 2 4

x1 + x2 + 2x3 = 2 A= 1 1 2

4x1 + 3x2 + 2x3 = 3 4 3 2

Usando o processo de Gauss, sem estratégia de pivoteamento parcial, para triangular A, temos:

Etapa1:⁽⁰⁾

A11(0) = 3 Multiplicadores: m21 = a21

(0) = 1 e m31 = a31(0) = 4 e

a11(o) 3 a11

(o) 3

Então,

L1 ⃪ L1 3 2 4

L2 ⃪ L2 -.m21L1 e A(1) = 0 1/3 2/3

L3 ⃪L3.- m31L1 0 1/3 -10/3

Uma vez que os elementos a21(1) e a31

(1) são nulos, podemos guardar os os multiplicadores nestas posições, então:

3 2 4

A(1) = 1/3 1/3 2/3

4/3 1/3 -10/3

Etapa 2:

Pivô: a22(1) = 1/3 Multiplicadores: m21 = a32

(1) = 1/3 = 1 a22

(1) 1/3 Teremos:

L1 ⃪ L1 3 2 4

L2 ⃪ L2. e A(1) = 1/3 1/3 2/3

L3 ⃪L3 -.m32L2 4/3 1 - 4

Os fatores L e U são:

1 0 0 3 2 4

L = 1/3 1 0 e U= 0 1/3 2/3

4/3 1 1 0 0 -4

Resolvendo L(Ux) = b:

i) Ly =b

y1 = 1

1/3y1 + y2 = 2

4/3y1 + y2 + y3 = 3

Y = ( 1 5/3 0)T

ii) Ux = y:

3x1 + 2x2 + 4x3 = 1

Ux = y → 1/3x2 +2/3x3 = 5/3

-4x3 = 0

X = (-3 5 0)T.

2.3.9 Métodos Iterativos para a Solução de Sistemas Lineares

Seja os Sistema Linear A x=b onde:

A matriz de coeficientes n×n

x vetor de variáveis n×1

b vetor independente (constantes) n×1

Idéia Geral dos Métodos Iterativos

Converter o sistema de equações A x=b em um processo iterativo

x=C x+g=ϕ ( x ) , onde:

C matriz com dimensões n×n

g vetor com dimensões n×1

ϕ ( x ) função de iteração matricial

Esquema Iterativo Proposto

Partindo de uma vetor aproximação inicial x(o )

, constrói-se uma seqüência iterativa de vetores:

x(1)=C x(o )+g=ϕ ( x(o ))

x(2)=C x(1 )+g=ϕ (x(1 ))

⋮

x(k )=C x(k−1)+g=ϕ( x(k−1))

Forma Geral

x(k+1)=ϕ ( x(k ))

Observação

Se a sequência de aproximação x(o )

, x(1)

, x(2)

, ......, x(k )

é tal que limk→∞

x(k )=α ⇒ α=C α+g, então α é a solução do sistema A x=b .

Teste de ParadaComo em todos os processos iterativos, necessita-se de um critério para a

parada do processo.

a) Máximo desvio absoluto:

Δ(k )=maxi=1, n

|x i(k )−x i

(k−1)|

b) Máximo desvio relativo:

ΔR(k )= Δ(k )

maxi=1, n

|x i(k )|

Desta forma, dada uma precisão ε , o vetor x(k )

será escolhido como solução

aproximada da solução exata, se Δ(k )<ε , ou dependendo da escolha, ΔR

(k )<ε .

3.5.1 Método Iterativo de Gauss-Jacobi

Considere o sistema linear:

{a11 x1+a12 x2+a13 x3+. .. .. . .. .. . .a1n xn=b1 ¿ {a21 x1+a22 x2+a23 x3+.. . .. .. . .. .. a2n xn=b2¿ {. .. .. . .. .. . .. .. . .. .. . .. .. . .. .. .. . .. .. . .. .. . .. .. . .. .. . .. .. . .. .. .. . .¿ {.. . .. .. . .. .. . .. .. . .. .. . .. .. .. . .. .. . .. .. . .. .. . .. .. . .. .. . .. .. .. . .. .. ¿ ¿¿¿

Supondo a ii≠0 , i=1,2 ,. . ., n , isola-se o vetor x mediante a separação

pela diagonal da matriz de coeficientes.

x1(k+1)=

1a11

(b1−a12 x2( k )−a13 x3

( k )−. .. . .. .. . .−a1n xn(k ))

x2(k+1)=

1a22

(b2−a21 x1(k )−a23 x3

(k )−. .. .. .. . ..−a2 n xn(k ))

⋮

xn(k+1)=1

ann

( bn−an1 x1(k )−an 2 x2

( k )−.. . .. .. .. .−ann−1 xn−1(k ) )

Assim, tem-se o sistema iterativo x=C x+g , onde:

C=[ 0 −a12 /a11 −a13/a11 ⋯ −a1n /a11

−a21/a22 0 −a23/a22 ⋯ −a2 n /a22

⋮ ⋮ ⋮ ⋮ ⋮−an 1/ann −an 2/ann −an 3/ann ⋯ 0 ]

g=[b1/a11

b2 /a22

⋮bn /ann

]

Dado uma aproximação inicial x(o )

, o Método de Gauss-Jacobi consiste em

obter uma seqüência x(o )

, x(1)

, x(2)

, ......, x(k )

, por meio da relação recursiva:

x(k+1)=C x(k )+g

Observe que o processo iterativo utiliza somente estimativas da iteração anterior.

Exemplo: Resolver o sistema de equações lineares, pelo Método de Gauss-Jacobi com

solução inicial x(o )= [0,7 −1,6 0,6 ]T e tolerância ε≤0 , 05 .

10 x1+2 x2+x3=7x1+5 x2+x3=−82 x1+3 x2+10 x3=6

Separando-se os elementos diagonais, tem-se:

x1(k+1)=

110(7−2 x2

(k )−x3(k ))

x2(k+1)=1

5(−8−x1

(k )−x3(k ))

x3(k+1)=

110(6−2 x1

(k )−3 x2(k ))

C=[ 0−2

10−1

10−1

50

−15

−210

−310

0 ] g=[7

10−8

56/10

]Solução para k=0 x

(1)=C x(0 )+g⇒ x(1)

x1=[0−2

10−1

10−1

50

−15

−210

−310

0 ] [0,7−1,60,6 ]+[710

−85

6 /10]

x1=[0 ,96−1 ,860 ,94 ]

Cálculo de ΔR(1)

:

|x1(1 )−x1

(0 )|=|0,7−0 ,96|=0 , 26

|x2(1 )−x2

(0 )|=|−1 ,86−1,6|=0 , 26

|x3(1 )−x3

(0 )|=|0,6−0 ,94|=0 , 34

ΔR(1)=0 , 34

maxi=1,3

|x i(1 )|=0 , 34

1 ,86=0 , 1828

Para k=1:

x(2)=[ 0 , 978−1 , 980 , 966 ] ⇒ ΔR

(2)=0 ,121, 98

=0 ,0606>ε

Para k=2:

x(2)=[ 0 , 9994−1 , 98880 , 99984 ] ⇒ ΔR

(3)=0 , 03241 , 9888

=0 , 0163<ε

x̄=[ 0 ,9994−1 , 9888

0 , 9984 ] é solução com erro menor que 0,05.

Condições Suficientes para a Convergência do Método de Gauss-Jacobi

Teorema

Seja o sistema linear A x=b e seja:

¿ j≠k ¿¿¿n¿¿ ¿|akk|

¿¿

Se α=max

k=1, nαk<1

, então o método G-J gera uma seqüência {x(k )} convergente para a

solução do sistema dado, independentemente da escolha da aproximação inicial x(0 )

.

Observe que esta é uma condição suficiente, se for satisfeita o método converge, entretanto se não for satisfeita nada se pode afirmar.

Exemplo 1:

Seja a matriz do exemplo dado anteriormente:

A=[10 2 11 5 12 3 10 ]

α1=(2+1)10

=0,3<1

α2=(1+1)

5=0,4<1

α3=(2+3 )10

=0,2<1

Tem-se a convergência garantida para qualquer vetor inicial.

Exemplo 2:

Seja o sistema de equações lineares:

x1+x2=3x1−3 x2=−3

α 1=11=1

α 2=13

As condições de convergência do teorema não são satisfeitas, entretanto o Método de Gauss-Jacobi gera uma seqüência convergente para a solução exata

x=[3 23

2 ]T

. Se as condições de suficiência não são satisfeitas, não significa que o método não possa convergir.

Exemplo 3:

Considere o sistema linear:

x1+3 x2+x3=−25 x1+2 x2+2 x3=30 x1+6 x2+8 x3=−6

A=[1 3 15 2 20 6 8 ]

α 1=(3+1)

1=4>1

α2=(5+2)

2=3,5>1

α 3=(0+6 )

8=0 ,75<1

As condições do teorema não são satisfeitas. Uma solução possível é permutar as equações. Seja no exemplo permutar a primeira equação com a segunda equação:

A=[5 2 21 3 10 6 8 ]

α 1=(2+2 )

1=0,8<1

α 2=(1+1 )

3=0 , 66<1

α3=(0+6 )

8=0 , 75<1

As condições passam a ser satisfeitas e a convergência é garantida para qualquer vetor inicial. Este tipo de procedimento nem sempre é possível.

Fórmula Matricial do Método Gauss-Jacobi

Decompõe-se a matriz de coeficientes A em:A=L+D+U

Onde:L – Matriz Triangular InferiorD – Matriz DiagonalU – Matriz Triangular Superior

L=[0 0 0 ⋯ 0

a21 0 0 ⋯ 0a31 a32 0 ⋯ 0

⋮ ⋮ ⋮ ⋱ ⋮an 1 an2 ⋯ ann−1 0

] D=[d11 0 0 ⋯ 00 d22 0 ⋯ 0

0 0 d33 ⋯ 0

⋮ ⋮ ⋮ ⋱ ⋮0 0 ⋯ 0 dnn

] U=[0 a12 a13 ⋯ a1 n

0 0 a23 ⋯ a2 n

0 0 0 ⋯ ⋮⋮ ⋮ ⋮ ⋱ an−1n

0 0 ⋯ 0 0]

(L+D+U )x=bL x+D x+U x=bD x=b−(L+U ) xD x(k+1 )=b−(L+U ) x(k )

x(k+1)=D−1 b−D−1(L+U )x(k )

3.5.2 Método Iterativo de Gauss-Seidel

Assim como no Método de Gauss-Jacobi o sistema linear A x=b é escrito na forma equivalente:

x=C x+gComo no Método Gauss-Jacobi, é realizada uma separação diagonal, e o

processo iterativo de atualização é seqüencial, componente por componente. A diferença é que, no momento de realizar-se a atualização das componentes do vetor numa determinada iteração, a formulação utiliza as componentes da iteração já atualizadas na iteração atual, com as restantes não atualizadas da iteração anterior. Por

exemplo, ao se calcular a componente x j(k+1)

da iteração (k+1), utiliza-se no cálculo as

componentes já atualizadas x1(k+1) , x2

(k+1 ) , . . . , x j−1(k+1 )

com as componentes ainda não

atualizadas da iteração anterior x j+1(k ) , x j+2

(k ) , . .. , xn(k )

.

x1(k+1)=

1a11

(b1−a12 x2(k )−a13 x3

(k )−a14 x 4(k )−.. . .. .. . ..−a1 n xn

(k ))

x2(k+1)=

1a22

(b2−a21 x1(k+1 )−a23 x3

(k )−a24 x4(k ). . .. .. . .. .−a2n xn

(k ))

x3(k+1)=1

a22

(b3−a31 x1(k+1 )−a32 x2

(k+1 )−a34 x 4(k )−.. .. . .. .. .−a2 n xn

(k ))

⋮

xn(k+1)=1

ann

( bn−an1 x1(k+1)−an2 x2

(k+1)−an3 x3(k+1 ) .. .. . .. .. .−ann−1 xn−1

(k+1))

Exemplo: Resolver o sistema linear utilizando o Método Iterativo de Gauss-Seidel, com

xo=[0,0,0 ]T e tolerância ε≤5×10−2.

5 x1+x2+x3=53 x1+4 x2+x3=63 x1+3 x2+6 x3=0O processo iterativo é dado por:

x1(k+1)=

15(5−x2

(k )−x3(k ))

x2(k+1)=1

4(6−3 x1

(k+1)−x3(k ))

x3(k+1)=

16(0−3 x1

(k+1)−3 x2(k+1))

Para k=0 e xo=[0,0,0 ]T

x1=[ 10 , 75−0 , 875]

Cálculo de ΔR(1)

:

|x1(1 )−x1

(0 )|=|1,0−0|=1,0

|x2(1 )−x2

(0 )|=|0 , 75−0,0|=0 ,75

|x3(1 )−x3

(0 )|=|−0 ,875−0|=0 ,875

ΔR(1)=

maxi=1,3

Δi(1)

maxi=1,3

|x i(1 )|=1,0

1,0=1,0>ε

Para k=1 e x1=[1,0 , 0 ,75 , −0 ,875 ]T :

x(2)=[1 , 0250 ,95−0 ,9875]

|x1(2 )−x1

(1 )|=|1,025−1,0|=0 , 025

|x2(2 )−x2

(1 )|=|0 , 95−0 ,75|=0 , 20

|x3(2 )−x3

(1 )|=|−0 ,9875−(−0 ,875 )|=0 ,1125

ΔR(2)=

maxi=1,3

Δi(2 )

maxi=1,3

|x i(2)|=0 , 20

1 ,025=0 , 1957>ε

Para k=2 e x2=[ 1,0025 , 0 , 95 , −0 , 9875 ]T :

x(3 )=[1, 00750 ,9912−0 ,9993]

|x1(3 )−x1

(2 )|=|1 ,075−1 , 025|=0 ,0175

|x2(3 )−x2

(2 )|=|0 , 9912−0 , 95|=0 ,0412

|x3(3 )−x3

(2 )|=|−0 ,9993−(−0 , 9875)|=0 , 0118

ΔR(3 )=

maxi=1,3

Δi(3 )

maxi=1,3

|xi(3)|=0 , 0412

1 , 0075=0 ,0408<ε

x̄=[ 1 ,00750 ,9912−0 , 9993] é solução com erro menor que 0,05.

Fórmula Matricial do Método Gauss-Seidel

Decompõe-se a matriz de coeficientes A em:A=L+D+U

Onde:L – Matriz Triangular InferiorD – Matriz DiagonalU – Matriz Triangular Superior

L=[0 0 0 ⋯ 0

a21 0 0 ⋯ 0a31 a32 0 ⋯ 0

⋮ ⋮ ⋮ ⋱ ⋮an 1 an2 ⋯ ann−1 0

] D=[d11 0 0 ⋯ 00 d22 0 ⋯ 0

0 0 d33 ⋯ 0

⋮ ⋮ ⋮ ⋱ ⋮0 0 ⋯ 0 dnn

] U=[0 a12 a13 ⋯ a1 n

0 0 a23 ⋯ a2 n

0 0 0 ⋯ ⋮⋮ ⋮ ⋮ ⋱ an−1n

0 0 ⋯ 0 0]

(L+D+U )x=bL x+D x+U x=bD x=b−(L+U ) xx=D−1b−D−1 L x−D−1U x

x(k+1)=D−1b−D−1 L x(k+1)−D−1U x(k )

Interpretação Geométrica do Caso 2×2

Considere o Sistema Linear:x1+x2=3x1−3 x2=−3O esquema iterativo utilizando o Método de Gauss-Seidel é dado por:

x1(k+1)=

11(3−x2

(k ))

x2(k+1)=1

3(3+3 x1

(k+1 ))

Para k=0 e x(0 )=[0 0 ]T :

x1(1)=3

x2(1)=2

Para k=1 e x(1)=[3 2 ]T :

x1(2)=1

x2(2)=4

3

Para k=2 e x(2)=[1 4

3]T

:

x1(3 )=5

3

x2(3 )=14

9

A solução exata é dada por: x=[3

23

2]T

.Esse processo iterativo até à convergência pode ser interpretado geométricamente num

grafico com a componente x1 na abscissa e a componente x2 na ordenada.

-5 -4 -3 -2 -1 0 1 2 3 4 5-2

-1

0

1

2

3

4

5

6

7

8x2

x1 (0,0) (0,3)

(3,2) (1,2)

(1,4/3) (4/3,5/3)

x1+x2=3

x1-3x2=-3

Observação 1: Verifica-se pelo gráfico acima que a seqüência x(o )

, x(1)

, x(2)

, ......, x(k )

está convergindo para a solução exata x=[3

23

2]T

.Observação 2: A sequência gerado pelo Método de Gauss-Seidel depende fortemente da posição das equações. Esta observação pode ser melhor entendida modificando a ordem das equações do exemplo anterior.

Considere o mesmo Sistema Linear anterior, porém modificando a ordem das equações:x1−3 x2=−3x1+x2=3O esquema iterativo utilizando o Método de Gauss-Seidel é dado por:x1(k+1)=(−3+3 x2

(k ))

x2(k+1)=(3−x1

(k+1))

Para k=0 e x(0 )=[0 0 ]T :

x1(1)=−3

x2(1)=−6

Para k=1 e x(1)=[−3 −6 ]T :

x1(1)=15

x2(1)=−12

Estudo da Convergência do Método de Gauss-Seidel

Existem dois critérios de suficiência para a convergência do Método de Gauss-Seidel. O critério de linhas e o critério de Sassenfeld. O critério de linhas é o mesmo da Método de Gauss-Jacobi.

-5 -4 -3 -2 -1 0 1 2 3 4 5-2

-1

0

1

2

3

4

5

6

7

8

x (0,0)

(0,-3)

(-3,6) ......(6,15)

x2

x1

x1+x2=3

x1-3x2= -3

Critério de Linhas

Seja o sistema linear A x=b , com A dimensão n×n e seja:

¿ j≠k ¿¿¿n¿¿ ¿|akk|

¿¿

Se α=max

k=1, nαk<1

, então o método Gauss-Seidel gera uma seqüência {x(k )} convergente para a solução do sistema dado, independentemente da escolha da

aproximação inicial x(0 )

.

A matriz que satisfizer o critério de linhas é chamada de diagonal dominante estrita.

Critério de Sassenfeld

Seja o sistema linear A x=b , com A dimensão n×n e seja:

β1=|a12|+|a13|+|a14|+ .. .. .. .+|a1n|

|a11|

e para j=2,3 ,. . .. .. . .. .. . .n :

β j=|a j 1|β1+|a j 2|β2+. .. .. . .. .. .. . .+|a jj−1|β j−1+|a jj+1|+.. . .. .+.|a j1|

|a jj|

Define-se β=max

j=1 , nβ j

.

Se β<1 , então o Método de Gauss-Seidel gera uma sequência convergente para a solução do sistema, qualquer que seja o vetor inicial. Além disso, quanto menor for o

valor de β mais rápida é a convergência.

Exemplo: Verificar as condições de convergência do Método de Gauss-Seidel no sistema abaixo:

2 x1+x2+3 x3=3−x2+x3=1

x1 + 3x3=3

a) Critério de Linhas

α 1=1+3

2=2>1

não satisfaz.b) Critério de Sassenfeld

β1=1+3

2=2>1

não satisfaz.

Como a convergência do Método de Gauss-Seidel é fortemente dependente da posição das equações, pode-se trocar a posição das equações.Tentativa 1: Troca-se a primeira equação pela terceira equação.x1 + 3x3=3−x2+x3=1

2 x1+x2+3 x3=3a) Critério de Linhas

α 1=0+3

1=3>1

não satisfaz.b) Critério de Sassenfeld

β1=0+3

1=3>1

não satisfaz.Tentativa 2: Troca-se a primeira coluna pela terceira coluna na equação anterior.3 x3+0 x2+x1=3x3−x2+0 x1=13 x3+x2+2 x1=3a) Critério de Linhas

α 1=13=0 . 33<1

satisfaz.

α 2=11=1

não satisfaz.

b) Critério de Sassenfeld

β1=13=0 . 33<1

satisfaz.

β2=1⋅1

31=0 .33<1

satisfaz.

β3=3⋅1

3+1⋅1

32

=46<1

satisfaz.

Com a última modificação o sistema passa a ser convergente para qualquer vetor inicial. Modificações desse tipo são puramente acadêmicas e são difíceis de serem realizadas em sistemas reais. Principalmente pelas dimensões dos problemas, resultando num grande esforço computacional, e das incertezas quanto a sua eficiência.

Exemplo : Verifique a convergência do sistema abaixo pelo critério de linhas e Sassenfeld.

2 x1+x2=3x1 + x2=4

Critério de Linhas

α 1=12

α 2=1 Não satisfaz

Critério de Sassenfeld

β1=12

β2=1×1

21

=12 Satisfaz

Exemplo : Verifique a convergência do sistema abaixo pelo critério de linhas e Sassenfeld.

x1+0,5 x2−0,1x3+0,1x 4=0,20,2 x1+x2−0,2x3−0,1 x4=−2,6−0,1 x1−0,2 x2+x3+0,2x4=1,0

0,1 x1+0,3 x2+0,2 x3+x4=−2,5

3.5.3 Método da Sobrerelaxação Sucessiva

x i(k+1)=x i

(k )+w ( x̂i(k+1)− xi

( k ))

x̂ i(k+1)= 1

aii

(bi−∑j<i

a ij x j(k+1 )−∑

j>i

aij x j(k ))

1<w<2 - Utilizado para aumentar a velocidade de convergência.0<w<1 - Utilizado em sistemas com dificuldades de convergência pelo Método de gauss Seidel.

Condições Necessária e Suficiente para Convergência do Método de Gauss-Jacobi e Gauss-Seidel

Teorema

Seja b∈ℜn

e A=(M−N ) ∈ℜné não singular. Se M é não-singular e o

raio espectral de M−1 N satisfaçaρ(M−1 N )<1 , então o processo iterativo definido por

M x(k+1)=N x(k )+b converge para x=A−1 b para qualquer vetor x(0 )

.

ρ(M−1 N )=max {|λ|: λ∈ λ(M−1 N )}

Para o Método de Gauss-Jacobi:(L+D+U )x=bL x+D x+U x=bD x=b−(L+U ) xD x(k+1 )=b−(L+U ) x(k )

M=DN=−(L+U )

Para o Método de Gauss-Seildel:(L+D+U )x=bL x+D x+U x=bL x+D x=b−U x(L+D) x(k+1)=b−U x(k )

M=(L+D )N=−U

Comparação dos Métodos de Solução de Sistemas Lineares A x=b

Métodos Diretos:

1. Processos finitos (convergência para qualquer sistema não-singular);2. Apresentam problemas com erros de arredondamento;3. Utiliza-se técnicas de pivoteamento para amenizar os problemas de arredondamento;4. O processo de triangularização destrói a esparsidade da matriz de coeficientes.

Técnicas de Esparsidade são utilizadas para amenizar o enchimento da matriz.

5. Redução de esforço computacional é conseguido em soluções para vetores independentes adicionais com a matriz de coeficientes mantida constante, com a utilização da fatoração LU;

6. Para matrizes cheias a solução requer n3

operações sem considerar o pivoteamento.

Métodos Iterativos:

1. Provavelmente mais eficientes para sistemas de grande porte, principalmente com a utilização de computação de alto desempenho (vetorial e paralela);

2. Tem convergência assegurada apenas sob certas condições;3. Conserva a esparsidade da matriz de coeficientes;

4. Os métodos de G-J e G-S requerem 2 n2operações por iterações;

5. Poucas vantagens adicionais são conseguidas em soluções para vetores independentes adicionais com a matriz de coeficientes mantida constante, como no caso da fatoração LU;

6. Carregam menos erros de arredondamento no processo, tendo em vista que a convergência uma vez assegurada independe da aproximação inicial. Somente os erros da última iteração afetam a solução.

2.4 Número de condicionamento de uma matriz Os elementos da matriz de coeficientes e do vetor independente de um sistema de equações lineares são na grande maioria das aplicações inexatos. Esta falta de exatidão pode ser originada porque os dados são resultantes de experimentos ou computados através de operações que carregam erros de arredondamento, ou mesmo do próprio armazenamento dos elementos em uma aritmética finita. A questão é quando a perturbação introduzida em elementos do sistema podem alterar a resposta. O algoritmo de eliminação de Gauss com pivoteamento pode ser considerado numericamente estável, o que pode-se assegurar que, para um sistema bem comportado, produzirá pequenos resíduos, mesmo para pequenas perturbações nos elementos do sistema. Portanto, alterações na resposta do sistema está associada ao comportamento do sistema. Este comportamento é medido pelo número de condição (condicionamento) da matriz.

Para entender o número de condicionamento de uma matriz é preciso relembrar o conceito de norma de vetores e matrizes.

Norma de vetores

A norma de vetores em Rn

é uma função || . || := Rn→ R que satisfaz

1. x≠0→‖x‖>0

2. ‖α x‖=|α| ‖x‖

3. ‖x+ y‖≤‖x‖+‖y‖

Embora existem uma infinidade de normas, as mais utilizadas na prática são:

Norma 1: ‖x‖1=∑

i

|x i|

Norma 2: ‖x‖ 2=√∑

i

x i2

Norma ∞ : ‖x‖∞=maxi|x i|

Normas de matrizes A norma de uma matriz é definida de forma análogo a norma de vetores. A norma de

uma matriz em Rm×n

é uma função || . || : Rmxn→R que satisfaz:

1. Α≠0⇒‖Α‖>0

2. ‖α Α‖=|α| ‖Α‖

3. ‖Α+B‖≤‖Α‖+‖B‖

Propriedade Adicional

‖ΑB‖≤‖Α‖‖B‖Sempre que o produto A .B é definido.

OBS: Uma norma de vetor || . ||v é consistente com uma norma de matriz || . ||M se :

‖Α x‖≤‖Α‖‖x‖

As principais normas utilizadas são:

‖Α‖1=Maxj∑i=i

n

‖a ij‖ (máxima Soma absoluta de colunas de A )

‖Α‖2=[máx auto valor de ΑT Α ]1

2= max valor singular de Α

‖Α‖∞=Maxi∑j=1

n

‖a ij‖ (máxima soma absoluto de linhas de Α )

Condição de uma matriz não-singular

Seja o sistema linear:Α x=b

Suponha que o vetor independente b sja alterado para b+δ b e a matriz A permanece inalterada. A solução exata do problema alterado é dado por:

Α (x+δ xb )=b+δ b

Como Α x=b , pode-se obter um limite para δ xb ;

‖δ xb‖≤‖Α−1‖|δ b|

A relação Α x=b implica em :

‖b‖≤‖Α‖‖x‖Multiplicando as duas relações anteriores chega-se a relação:

‖δ xb‖‖x‖

≤‖Α−1‖‖Α‖‖δ b‖‖b‖

Se perturbarmos a matriz A enquanto b é mantido fixo tem-se a solução: ( Α+δΑ ) (x+δxΑ )=b

por um processo similar, chega-se a expressão para o erro relativo:‖δx Α‖‖x+δx Α‖

≤‖Α−1‖‖Α‖‖δΑ‖‖Α‖

A quantidade ‖Α−1‖‖Α‖ reflete a máxima mudança relativa possível na solução exata para um sistema linear com dados perturbados, portanto:

cond ( Α )=‖Α−1‖‖Α‖Das relações anteriores pode-se chegar a:

cond ( Α )≥

‖δxb‖||x||

‖δb‖||b||

e cond ( Α )≥

‖δxΑ‖‖x+δx Α‖‖δΑ‖‖Α‖

Estas relações mostram que, se a mudança relativa é muito grande para uma pequena

perturbação em A ou b , nós sabemos que A é mal condicionada.

Propriedades:

1. A norma da matriz identidade, para qualquer norma, vale 1.

2. Como I=A−1 A e ‖A−1 A‖≤‖A−1‖‖A‖ conclui-se que cond ( A )≥1 .

3. Se a matriz A for multiplicada por um escalar α , então cond ( αA )=cond (A ) .

4. Se D for uma matriz diagonal, então cond (D )=

max‖d ii‖min‖d ii‖ .

5. Se A for não-singular e simétrica cond ( A )=

λmax

λmin . Se A for não-simétrica

cond ( A )≥λmax

λmin .

6. Se A for não-singular , Cond (A )=

σ1

σ n , onde σ 1é o mínimo valor singular e σ n é o máximo valor singular.

A=UΣV T onde

Σ=diag {σ i}i=1 , . .. , n .

Observação 1: Ocond ( A ) pode ser uma medida mais eficaz para a verificação da singularidade de uma matriz, que o determinante. Como exemplo, seja uma matriz

diagonal de ordem 100×100 , com todos elementos igual a 0,1. O determinante da

matriz é 10−100, que é um número pequeno e pode ser arredondado para zero em

computadores digitais. Já a número de condição da matriz é igual a 1.

Observação 2: O pré-condicionamento de uma matriz está associado na redução do raio

espectral da mesma.. ς= {autovalor de maior valor absoluto }

Α x=bΡΑ x=Ρ b

A matriz A é pré-multiplicada por uma matriz P de pré-condicionamento, com o objetivo de que a nova matriz PA tenha uma redução do raio espectral o quê implica na redução do condicionamento.

O mal condicionamento de uma matriz está associado à proximidade da singularidade da matriz. Exemplo 2x2.

{ax1+bx2=cdx1+cx2=f

Matriz não-singular

Matriz singular

Matriz próxima da singularidade

ExemploΑ x=bA matriz A é armazenada em uma máquina com unidade de arrendondamento u .

~a ij=aij (1+εij ) [ εij ]≤u~Α=Α+Ε onde e ij=aij εijPara qualquer norma:

‖Ε‖≤‖Α‖USe o algoritmo de solução for numericamente estável:‖~x−x‖‖~x‖

≤ cond ( Α )u

Condicionamento de Matrizes de Redes Elétricas

Análise Nodal A equação matricial de desempenho de uma rede linear na formulação nodal é dada por: I=YEonde I – é o vetor das injeções nodais de corrente E – vetor das tensões nodais em relação ao nó de referência. Y – matriz de admitância nodalEm sistemas de potência geralmente o nó de referência é o nó terra, sendo os demais nós as barras do sistema: IBarra = YBarra ZBarra

Na sua forma inversa tem-se: EBarra = YBarra

-1 IBarra

YBarra → esparsa

ZBarra →Cheia

YBarra-1 →ZBarra

Na ausência de acoplamento mútuo, a matriz YBarra é montada com grande facilidade.

Elementos diagonais Υ ii=∑

K

y ik

K barra vizinha a i e y ik admitância da linha i-k

Υ iK = { − y ik ∀ K vizinho a i0 ∀K não vizinho a i

Exemplo

1 2

I1 1 - j2 E

-j4

1 2 – j4

4 3

I4

I3

J1/2 -J

110

[ Ι 1

−J 4 E

−I 3

−I 4]=[ 2− j 2 −1+ j 2 0 −1−1+ j2 3− j10 −2+ j 4 0

0 −2+ j 4 3− j7,1 −1+ j3−1 0 −1+ j 3 2− j 2,5

] [Ε1

Ε2

Ε3

Ε4]

Condicionamento da matriz YBarra

Sistema sem ligação à terra Yb

11

Ya Yc

YBarra = [Ya+Yb −Ya −Yb−Ya Ya+Yc −Yc−Ya −Yc Yb+Yc ]

Matriz é singular (colunas combinação linear)Mau condicionamento: admitância fraca entre a rede e o nó de referênciaBom condicionamento: forte conexão com o nó de referência.

~

-J4 E

-J4

EQUIVALENTE NORTON

I=VR

I=YV

13

Possibilidades de Mau Condicionamento1. Conexão fraca com o nó de referência2. Junção entre partes com admitâncias muito grandes e pequenas (perda dominância diagonal, aumento do raio espectral)3. Capacitância em série ou em derivação do sistema, enfrequecendo o elemento diagonal.

Melhoria do Condicionamento adotando como referência uma barra do sistema. Se a tensão de uma barra for considerada conhecida, pode-se reduzir o n° de equações e variáveis de uma unidade. I = Y E referência ao nó terraΙ 1=Υ 11 Ε1 + Υ 12 Ε2 +.. . .. .. . .Υ 1 n Εn

Ι 2=Υ 21 Ε1 + Υ 22 Ε2 +. . .. .. . .. Υ 2n Εn

Ι n=Υ n1 Ε1+Υ n 2 Ε2+. . .. .. . .. ..Υ nn Εn

Com Ε1 de referênciaΙ 1−Υ 11 Ε1=Υ 12 Ε2+. . .. .. . .+Υ 1n Εn

Ι 2−Υ 21 Ε1=Υ 22 Ε2+. .. .. . .. ..+Υ nn Εn

Ι n−Υ n1 Ε1=Υ n2 Ε2+.. . .. .. . .. .+Υ nn Εn

Pode-se eliminar uma das equações, já que há apenas n-1 incognitas.

Ι̂−C=Υ̂ Ε̂Υ̂ obtida de Y eliminando-se a linha e coluna correspondentes à nova barra de

referência. Para obter-se bom condicionamento de Υ̂ escolhe-se como referência uma barra com forte conexão à rede restante.

2.5 Correção IterativaPara sistemas de equações lineares que não se conhece a priori se é bem

condicionado é importante verificar se a solução é suficientemente exata; e se sua exatidão não for suficiente, pode-se melhorá-la.

Seja o sistema Linear:Α x=b

A exatidão da solução desse sistema pode ser verificada pelo resíduo:

r(1 )=b−Α x (1 )

x (1 )

solução computada.

Se os elementos de r são pequenos, se comparados aos elementos de b , normalmente assume-se que a solução é suficientemente exata.

Se a solução não for suficientemente exata, pode-se repetir a solução em dupla precisão ( o que normalmente não acrescente grandes beneficios, principalmente se a matriz de coeficientes for mal condicionada).

Um procedimento que pode ser adotado para melhorar a exatidão é a correção iterativa, conforme a sequência:

r (1 )=b−Α x(1)

b=Α x

Α⋅( x−x (1 ) )=r (1 )

y(1 )=x−x (1 )

Chega-se a um novo sistema linear:

Α y(1 )=r(1)Resolvendo-se esse sistema, uma correção da variável é obtida:

x(2)=x(1)+ y(1)Esse processo pode ser repetido até que se chegue a convergência.

Observações:

a) A decomposição da matriz é realizada uma única vezb) Em razão do mal condicionamento do sistema, os arredondamentos de elementos da matriz de coeficientes A podem causar grandes erros na solução, é muito importante que a matriz de coeficientes seja construída em dupla precisão e armazenada em dupla precisão. É também necessário computar o resíduo em dupla precisão, de forma que erros significativos não sejam introduzidos na computação.c) É possível que o método não convirja, se a amplificação do erro for muito grande.

e (k )=

‖r (k )‖‖b‖

Observe que:

Uma saída deve ser prevista a implementação do processo no caso em que e(k+1 )>e (k )

Início

Entre com A e b em dupla precisão

2.6 Eliminação por Blocos Aproveitar estrutura esparsa de matrizes por blocos.

Fluxo de Potência ótimo (Mesma estrutura da matriz admitância, se for feito por blocos) Fluxo de Potência via Newton Raphson.

Copie A em precisão simples

Decomponha a Matriz A em fatores LU

Faça k=0 , x(0 )=0 (precisão dupla),

Solucione o Sistema A y(k )=r(k ) (precisão simples)

Corrija a Solução:

x(k+1)=x(k )+ y(k ) (precisão dupla)

Teste

e(k+1)≥e(k )

Saída com Solução Insatisfatória

e(k+1)≤tolerância

Saída com Solução Satisfatória

e(k+1)<e(k )

k=k+1

Α=(Α11 Α12 . . .. Α1 m

Α21 Α22 . . .. Α2m

Αm 1 Αm 2 . . .. Αmm

)Αii são blocos diagonais quadrados de ordem ni

n=n1+n2+. .. . .. .. . nm

A maneira mais simples é realizar através de computação recursiva:

Α=(Α11 Α1∗¿ ¿Α¿1 Α** )

multiplicadores: Α¿ 1 . Α11−1

Α=( I 0Α ¿1 Α11

−1 I ) ¿¿Explicitando.1. A primeira linha de blocos de U é a primeira linha de blocos de A.

2. A primeira coluna de blocos de L excluindo o bloco identidade diagonal é Α¿ 1 Α11−1

.

3. Os blocos restantes da decomposição são obtidos por Α**−Α¿1 Α11−1 Α1∗¿ ¿

A decomposição escalar e a decomposição por blocos normalmente não são iguais.

O número de operações é o mesmo da eliminação escalar: Ο(n3

3)Algoritmo

Α=(IL21 ILm1 Lm2 I ) (

U11 U12 U1 m

0 U22 U2 m

0 U mm)

nn=n (1 )+.. .. . .. n(m ) → n ( i ) dim ensão de cada bloco .ux=0Para K=1 a m−1

lx=ux+1ux=lx+n (K )−1if ( A [ lx :ux , lx :ux ] é singular ) erro

A [ux+1 :nn , lx :ux ]=A [ux+1 :nn , lx :ux ]∗A [ lx :ux , lx :ux ]−1

A [ux+1 :nn , ux+1: nn ]=A [ux+1 :nn , ux+1: nn ]−−A [ux+1:nn , lx :ux ]∗A [ lx :ux ,ux+1 :nn ]

Fim K