Vectores e Matrizes Aplicações à Engenharia Pedro Barahona DI/FCT/UNL Introdução aos Computadores e à Programação 2º Semestre 2007/2008 Tipo de Dados Primitivo: a Matriz • Ao contrário da generalidade das linguagens de programação, o Octave / MATLAB assume a matriz como o tipo básica de uma variável. >> V = 5 >> V(1) ans = 5 • Por exemplo, ao se fazer a atribuição de um valor a uma variável V, “simples”, o Octave está de facto a atribuir esse valor a todos os elementos de uma matriz com uma linha e uma coluna. >> V(1,1) ans = 5 • Para obtermos a dimensão de uma variável X, quer seja simples, vector ou matriz, podemos usar as funções prédefinidas rows(X) e columns(X). No caso de vectores, linha ou coluna, o número de colunas e linhas é retornado pela função lenght(V). >> X = 5 26 Março 2008 >> columns(X) ans = 1 >> M = [1 ; 2]; >> cols(M) ans = 2 Vectores e Matrizes - Aplicações à Engenharia 2 Matrizes • De facto, as funções aplicáveis a uma variável “simples” são sempre distribuídas por todos os elementos de uma matriz. Por exemplo, se a função logaritmo receber como parâmetro um vector ou matriz, retorna uma estrutura idêntica como resultado. >> V = [1,2, 4;8 16 32] V = 1 2 4 8 16 32 >> X = log(V)/log(2) X = 0 1 2 3 4 5 • Naturalmente, esta atribuição pode ser feita elemento a elemento. Assim, sendo V a anterior matriz, a instrução log(V)/log(2) é equivalente ao programa for i = 1:rows(V) for j = 1:columns(V) X(i,j) = log(V(i,j))/log(2) endfor endfor 26 Março 2008 Vectores e Matrizes - Aplicações à Engenharia 3 União e Seleccção de Matrizes • O facto de a estrutura matriz ser uma primitiva da linguagem Octave permite que a “formação” de matrizes seja feita não só a partir de elementos simples, como também o seja a partir de vectores ou matrizes. • Assim, a partir de dois vectores 1*3, A = [1,2,3] e B = [4,5,6], pode constituir-se % uma matriz de 2*3 >> C = [A;B] C = 1 2 3 4 5 6 • % um vector de 1*6 >> D = [A,B] C = 1 2 3 4 5 6 Similarmente, a selecção de elementos de uma matriz não tem de ser feita elemento a elemento, podendo ser seleccionadas submatrizes mais complexas . Por exemplo, >> V = [1,2,3; 4 5 6; 7 8 9; 3 6 9] V = S 26 Março 2008 1 4 7 3 2 5 8 6 3 6 9 9 T >> S = V(2:3, 1:2) S = 4 5 7 8 >> T = V(:,2:3) T = 2 3 5 6 8 9 6 9 Vectores e Matrizes - Aplicações à Engenharia 4 Selecção de Sub-Matrizes • De notar que esta selecção pode ser feita elemento a elemento, através de 2 ciclos encadeados, havendo necessidade de proceder a uma conversão de índices. Para a matriz anterior de 4*3, a instrução S = V(2:3, 1:2) é equivalente ao programa for i = 1:rows(V) for j = 1:columns(V) if i >= 2 & i <= 3 & j >= 1 & j <= 2 s(i-1,j) = v(i,j) endif endfor endfor • • De notar a utilização de expresssões booleanas compostas, obtidas pela aplicação de operadores booleanos (& - “e”, | - “ou” e ! – “não”) a relações de comparação simples (==, >=, <=, >, <, !=). 1 2 3 V De notar ainda a mudança de índices entre a matriz V e 4 5 6 a sua sub-matriz S: i(S) = i(V)+1; j(S) = j(V). 7 8 9 S 3 6 9 26 Março 2008 Vectores e Matrizes - Aplicações à Engenharia 5 Funções de Agregação sobre Vectores • Por vezes estamos interessados em conhecer propriedades do conjunto dos elementos de um vector, tais como os seus máximo, mínimo ou soma. • Estes valores agregados podem ser calculados por funções que utilizam ciclos, em que se vai alterando uma variável de acumulação. Por exemplo, para calcular o máximo ou a média de um vector, podem utilizar-se as seguintes funções function m = maximo(V); m = -inf; for j = 1:length(M) m = max(m, V(j)); endfor endfunction function s = soma(V); s = 0; for j = 1:length(M) s = s + V(j); endfor endfunction • Na realidade não é necessário definir estas funções, pois já existem predefinidas as funções max(V), min(V), sum(V), prod(V), mean(X). • Estas e outras funções agregadas podem ser obtidas indirectamente, como é o caso da média dos valores de um vector V: mean(X)= sum(V)/length(V). 26 Março 2008 Vectores e Matrizes - Aplicações à Engenharia 6 Soma de Vectores • Em geral, qualquer operação “” entre valores numéricos pode ser estendida a todos os valores corespondentes de vectores e matrizes antecedendo o operador por um ponto, i.e. através do operador “.”. • Desta forma a soma de dois vectores V1 e V2, que implementa a conhecida regra do “paralelograma”, pode ser especificada através da operação V = V1 .+ V2 • Por ser a matriz um tipo de dados primitivo, e por estar definida a soma de matrizes, a mesma operação pode ser, no caso da soma, ser especificada como V = V1 + V2 que faz o “overload” do operador "+”. V1 V >> V1 = [5,4] >> V2 = [7,-3] >> V = V1 + V2 ans = [12 1] V2 26 Março 2008 Vectores e Matrizes - Aplicações à Engenharia 7 Produto Interno de Vectores • Para além da soma de vectores, o produto interno de vectores é igualmente muito utilizado em Engenharia, por exemplo, para determinação do ângulo entre dois vectores e do módulo (valor absoluto, intensidade) de um vector. • Formalmente o produto interno entre dois vectores com o mesmo número de elementos, A e B, (denotado por A B) é a soma do produto de todos os elementos correspondentes. Pode pois ser especificado pela função prod_int abaixo (que pressupõe 2 vectores linha com o mesmo número de colunas) function p = prod_int(V1, V2); p = 0 for j = 1:columns(V1) p = p + V1(j)* V2(j); endfor endfunction • Como será de esperar, este produto interno pode ser obtido directamente através da composição das funções sum e o produto elemento a elemento, isto é X = prod_int(V1,V2) 26 Março 2008 é equivalente a X = sum( V1 .* V2) Vectores e Matrizes - Aplicações à Engenharia 8 Médias Ponderadas • É muito vulgar pretender-se obter a média de um conjunto de valores, ponderada pela sua “importância”. Por exemplo, a nota final é a média ponderada entre a nota prática e a nota escrita (exame), sendo os pesos relativos, 25% e 75%, respectivamente.. • Como é sabido, a média de k valores, v1 .. vk, com pesos relativos p1 .. pk define-se i n i n através da expressão m i 1 v(i ) * p(i ) / i 1 p(i ) • Dados os vectores V e P, contendo respectivamente os valores e os pesos, a média ponderada pode ser obtida através da função function m = media_ponderada(V, P); x = 0; s = 0; for j = 1:columns(V) x = x + V(j)* P(j); s = s + P(j) endfor m = x/s endfunction ou simplesmente através da expressão 26 Março 2008 m = sum(V .* P)/sum(P) Vectores e Matrizes - Aplicações à Engenharia 9 Módulo e Ângulo entre Vectores • Para além da soma de vectores, o produto interno de vectores é igualmente muito utilizado em Engenharia. • Por exemplo, o módulo de um vector A = [a1, a2, a3, ... , an], é definido como |A| = (a12 + a12+ ... an2)1/2 e portanto, a função vec_mod pode ser definida como function p = mod_vec(V); p = sqrt(prod_int(V,V)) endfunction • Para definir, na função ang_vec, o ângulo entre dois vectores A e B, basta notar que A B = | A | | B | cos(), e portanto = arccos(A B / (| A | | B | )), donde function a = ang_vec(V1,V2); % retorna ângulo em graus M = prod_int(V1,V2)) cos = M /(mod_vec(V1) * mod_vec(V2)) a = acos(cos)*180/pi endfunction 26 Março 2008 Vectores e Matrizes - Aplicações à Engenharia 10 Módulo e Ângulo entre Vectores • Exemplo: Um corpo é submetido a duas forças, F1 e F2. – F1 = 2 ex -3 ey+ 4 ez ; F2 = 1 ex + 2 ey - 2 ez Qual a intensidade (módulo) de cada uma das forças, qual o ângulo que essas forças fazem entre si, e qual a força resultante, bem como a sua intensidade (Módulo). F1 = [2 -3 4] , F2 = [1, 2, -2] >> M1 = mod_vec(F1,F1), M2 = mod_vec(F2) M1 = 5.3852 M2 = 3.0000 >> F F = = F1 + F2, 3 -1 M = mod_vec(F) 2 % vector [3 -1 2] M = 3.7417 >> alfa = ang_vec(F1,F2) % em graus alfa = 137.97 26 Março 2008 Vectores e Matrizes - Aplicações à Engenharia 11 Filtros • Tipicamente, as expressões booleanas são usadas para decidir condições de entrada e saída de execuções condicionais e/ou iterativas. • Na realidade, as expressões booleanas são avaliadas no conjunto {FALSE, TRUE} que no Octave coincide com o conjunto {0,1} de valores numéricos, o que permite trocar execuções condicionais por “filtros”. • Exemplo: Dado um valor numérico x, atribuir a y esse valor, caso seja maior que 5, ou 0, no caso contrário. Este problema pode ser resolvido através de (pelo menos) 2 maneiras distintas. Execução Condicional 26 Março 2008 if x > 5 y = x else y = 0 endif; Filtro y = x * (x > 5) Vectores e Matrizes - Aplicações à Engenharia 12 Filtros em Vectores • A filtragem de elementos pode ser feita naturalmente para todos os elementos de um vector, utilizando a extensão do operador de filtragem (geralmente a multiplicação, “*”) pelo correspondenete para todos os elementos do vector (i.e. “.*). • Exemplo: Determinar a soma de todos os elementos de um vector que sejam pares mas não múltiplos de 6. 1. Começamos por definir as funções booleanas par e múltiplo de 6. function p = div_2(x) p = (rem(x,2) == 0) endfunction; 2. Function p = div_6(x) p = (rem(x,6) == 0) endfunction; Podemos agora utilizá-las em expressões/filtros booleanos function s = sum_esp(X) p = sum(X .* (div_2(X) & !div_6(X)) endfunction; 26 Março 2008 Vectores e Matrizes - Aplicações à Engenharia 13 Matrizes • As operações descritas para vectores são na generalidade extensíveis a matrizes, com várias linhas e colunas. Por exemplo, – a matriz A pode ser “dobrada” através da operação A*2 – duas matrizes A e B com o mesmo número de linhas e colunas podem ser somadas quer através da operação A+B quer através de A.+B. • Existem no entanto algumas diferenças nas operações de agregação (max, min, sum, e prod) que não retornam os valores agregados de todos os elementos da matriz, mas os agregados, coluna a coluna. • Assim para se obter os valores agregados de toda a matriz deverá repetir-se a operação desejada, primeiro para agregação das colunas e depois para agregação da linha resultante. A = [ 1 2 3; 8 5 2]; >> sum(A) ans = 9 7 5 >> sum(sum(A)) ans = 21 26 Março 2008 Vectores e Matrizes - Aplicações à Engenharia 14 Multiplicação de Matrizes • Sendo a matriz o tipo de dados primitivo, o Octave implementa, como operação primitiva a multiplicação de matrizes. Dadas duas matrizes A e B, com dimensões m*k e k*n, o resultado é uma matriz C de dimensões m*n com elementos c(i, j ) l 1 a(i, l ) * b(l , j ) l k • Naturalmente, esta operação poderia ser igualmente obtida pela função mult_mat function C = mult_mat(A,B) for i = 1: rows(A) for j = 1: columns(B) C(i,j) = 0; for h = 1: columns(A) C(i,j) = C(i,j)+ A(i,h)*B(h,j); endfor endfor endfor endfunction 26 Março 2008 Vectores e Matrizes - Aplicações à Engenharia 15 Multiplicação de Matrizes • Exemplo: >> A = [ 1 2 3 4 ; 7 5 3 1 ; 0 2 4 6] A = 1 2 3 4 7 5 3 1 0 2 4 6 >> B = [ 1 2; 3 4; 5 6 ; 7 8] B = 1 2 3 4 5 6 7 8 >> C = C = 50 44 68 A*B 60 60 80 % 50 = 1*1 + 2*3 +3*5 + 4*7 >> D = (C == mult_mat(A,B)) D = 1 1 1 1 1 1 26 Março 2008 Vectores e Matrizes - Aplicações à Engenharia 16 Transposição de Matrizes • Em engenharia, a multiplicação de matrizes tem várias aplicações práticas. Para além da resolução de sistemas de equações, pode ser ainda utilizada na implementação dos produtos interno e externo de vectores (e ainda para rotação de vectores). • Para o produto interno, interessa definir inicialmente a operação de transposição de matrizes, que informalmente corresponde a trocar as linhas pelas colunas. Pode ser especificada formalmente através da função transp(A) definida abaixo. • Em Octave esta operação é primitiva, e é denotada pelo operador ‘, posfixo, ou seja A’ é equivalente a function B = transp(A) for i = 1: rows(A) for j = 1: columns(A) B(j,i) = A(i,j); endfor endfor endfunction 26 Março 2008 transp(A) A = 1 2 3 4 5 6 7 8 A’= 1 2 3 4 5 6 7 8 Vectores e Matrizes - Aplicações à Engenharia 17 Produto Interno • Para se obter o produto interno de dois vectores linha V1 e V2 com k elementos, basta notar que este produto pode ser obtido através da multiplicação da “matriz” V1 (1*k) com a “matriz transposta” de V2 (k*1). >> A = [ 1 3 5] A = 1 3 5 >> P = [ 5 2 1] P = 5 2 1 >> C = A * P’ C = 16 • Igualmente a média dos valores de um vector V, ponderada pelos pesos de outro vector P, pode ser obtida através da operação V * P’/sum(P) >> Mp = A * P’/sum(P) Mp = 2 26 Março 2008 Vectores e Matrizes - Aplicações à Engenharia 18 Produto Externo de 2 Vectores • Um outro produto de vectores (definido em espaços 3D) usado em Engenharia é o produto externo de dois vectores V1 e V2, V1 V2, definido como o vector V, com módulo igual ao produto dos módulos dos vectores V1 e V2 pelo seno do ângulo formado entre eles, com direcção perpendicular a ambos os vectores e com sentido definido pela regra do “saca-rolhas”. 1.8226 1.8226 = abs(1*2*sin(120*pi/180)) 2 120º 1 • Em particular este produto é usado para determinar a força exercida por uma carga eléctrica sujeita a uma força magnética, tal como acontece nos tubos de raios catódicos (CRT) usados nos “antigos” monitores de televisão. • Para obtermos este produto é conveniente a utilização de matrizes, da forma seguinte.. 26 Março 2008 Vectores e Matrizes - Aplicações à Engenharia 19 Produto Externo de 2 Vectores • Denotando por x, y e z os vectores unitários dos repectivos eixos, e de acordo com a definição, temos z – x x = y y = z z = 0 (pois x faz um ângulo de 0º com x – x y = z ; x z = - y; y z = x (regra do saca-rolhas) ; – y x = - z ; z x = y ; z y = - x (regra do saca-rolhas) x y • Sendo o produto externo distributivo em relação à soma, dados dois vectores A = ax x + ay y + az z e B = bx x + by y + bz z o seu produto externo é dado por • A B = (ax x + ay y + az z ) ( bx x + by y + bz z) = (ay bz - az by) x + (az bx - ax bz) y + (ax by - ay bx) z que pode ser obtido pela multiplicação do vector A pela matriz M abaixo indicada [ax ay az] 26 Março 2008 0 -bz by bz 0 -bx -by bx 0 = [aybz-azby azbx- axbz azbx-axbz] M Vectores e Matrizes - Aplicações à Engenharia 20 Produto Externo de 2 Vectores • O produto externo de dois vectores pode pois ser definido através da função function P = prod_ext(A,B) M = zeros(3,3); M(1,2) = -B(3); M(1,3) = B(2); M(2,1) = B(3); M(2,3) = -B(1); M(3,1) = -B(2); M(3,2) = B(1); P = A * M; endfunction • 0 -bz bz -by by 0 -bx bx 0 Por ser muito utilizado, o produto externo (cross-product em inglês) é disponibilizado em Octave pela função primitiva cross. Assim, >> A = [ 1 -3 5]; B = [ -1 0 4]; >> C = cross(A,B) C = -12 -9 -3 >> D = prod_ext(B,A) D = 12 9 3 26 Março 2008 % A B = - B A Vectores e Matrizes - Aplicações à Engenharia 21 Matrizes e Sistemas de Equações • Uma outra aplicação muito importante de matrizes em engenharia (e não só) é na resolução de sistemas de equações lineares. Um sistema de n equações lineares a n incógnitas, pode ser representado na forma matricial por a11 x1 + a12 x2 + ... + a1n xn = b1 a21 x1 + a22 x2 + ... + a2n xn = b2 AX=B ..... an1 x1 + an2 x2 + ... + ann xn = bn em que A é uma matriz n*n, e X e B são vectores coluna com n elementos, isto é ... a1n x1 b1 a22 ... a2n x2 b2 a11 a12 a21 ..... an1 26 Março 2008 an2 ... ann * ... xn = ... bn Vectores e Matrizes - Aplicações à Engenharia 22 Matrizes e Sistemas de Equações • Um tal sistema ficará resolvido se se colocar na forma 1 x1 + 0 x2 + ... + 0 xn = s1 0 x1 + 1 x2 + ... + 0 xn = s2 M*A*X = M*B ..... 0 x1 + 0 x2 + ... + 1 xn = sn I*X = S o que pode ser obtido através da multiplicação de uma matriz M em ambos os lados da equação. Como M * A = I, a matriz M é a matriz inversa da matriz A. m11 m12 ... m1n a11 a12 ... a1n 1 0 ... 0 m21 m22 ... m2n a21 a22 ... a2n 0 1 ... 0 ..... mn1 mn2 ... mnn • * ..... an1 an2 ... ann = ..... 0 0 ... 1 De notar que a multiplicação de matrizes não é comutativa, pelo que a multiplicação de B por M “à esquerda”, M * B, é diferente da multiplicação “à direita”, B * M, ou seja M*B ≠ B*M 26 Março 2008 Vectores e Matrizes - Aplicações à Engenharia 23 Matrizes e Sistemas de Equações • Existem várias formas de inverter uma matriz. Sendo a matriz um tipo de dados primitivo no Octave a operação de inversão de uma matriz A pode ser invocada da forma “standard”, A-1, ou alternativamente pela chamada da função pré-definida inv(A), ou seja A-1 • == inv(A). Assim, para obter as soluções de um sistema de equações A*X = B, basta notar que se A*X = B então é A-1*A*X = A-1*B, ou seja I*X = A-1*B, e finalmente X = A-1B. • Na álgebra “real”, a-1 * b = b/a, isto é a multiplicação pelo inverso corresponde a uma divisão. A divisão é igualmente uma operação primitiva do Octave. No entanto como as multiplicações “à direita” e “à esquerda” são diferentes, o Octave distingue a divisão “à esquerda” da divisão “à direita”, sendo representadas como A-1*B = A\B 26 Março 2008 e B*A-1 = B/A Vectores e Matrizes - Aplicações à Engenharia 24 Matrizes e Sistemas de Equações • Exemplo: Consideremos o sistema de 3 equações a 3 incógnitas 2 2 x1 + 4 x2 – x3 = 7 x1 - 2 x2 + x3 = 0 -3 x1 + 3 x2 – x3 = 2 • 1 -2 -3 x1 4 -1 1 4 -1 * x2 7 = x3 0 2 Para o resolver basta obter a matriz inversa A-1 e multiplicá-la à direita por B >> A = [2 4 -1; 1 -2 1 ; -3 3 -1] A = 2 4 -1 1 -2 1 -3 3 -1 >> M = A^-1 M = 0.20000 0.20000 0.20000 0.00000 -0.20000 0.50000 0.30000 2.00000 0.80000 >> B = [7 ; 0 ; 2] ; >> X = M*[7;0;2] X = 1 2 3 26 Março 2008 Vectores e Matrizes - Aplicações à Engenharia 25 Matrizes e Sistemas de Equações 2 2 x1 + 4 x2 – x3 = 7 x1 - 2 x2 + x3 = 0 -3 x1 + 3 x2 – x3 = 2 • 1 -2 -3 x1 4 -1 1 4 -1 * x2 7 = x3 0 2 Podemos ainda notar igualmente que >> A = [2 4 -1; 1 -2 1 ; -3 3 -1]; B = [7 ; 0 ; 2]; >> M = A ^-1; M = 1.00000 0.00000 0.00000 D = M * A 0.00000 0.00000 1.00000 0.00000 0.00000 1.00000 >> X = A\B X = 1 2 3 >> A \ B - M*B ans = 1.0e-16 * 4.44089 4.44089 4.44089 26 Março 2008 % erros de arredondamento! % (A \ B == M * B) = 0 0 0 Vectores e Matrizes - Aplicações à Engenharia 26