A SciPy Linalg

Veja esta apresentação em www.bordoni.info.

Transcript

1. A scipy.linalg Prof. Paulo R. G. Bordoni UFRJ
2. Aulas atrás, estudamos o básico da NumPy e da MatPlotLib. Iniciaremos agora nosso estudo da SciPy – Scientific Python.
3. E o que é, precisamente, Computação Científica Mestre? É o que você aprenderá ao longo deste curso, Loirinha. Vou começar pelo aspecto histórico!
4. http://history.siam.org/ A visão da Sociedade de Matemática Aplicada e Industrial – SIAM:
5. Procure por Scientific Computing na Web, Surfista. Farei a busca com o Google, Mestra.
6. Peguei o retorno do Google e separei o que achei interessante:
7. Comecei pela Wikipedia. Eis o que achei:
8. Nosso Mestre já citou Maple, Mathematica e MatLab. Como ele, a Wikipedia também indica Python e SciPy!
9. O 2º que separeifoi o site da INTEL de HPC. HPC é High Performance Computation
10. Separei algumas chamadas desse site que me interessaram – estão todas na mídia recente!
11. A SciPy está em 5º lugar na lista. E é reconhecida como uma biblioteca fundamental para computação científica! É com ela que trabalharemos Galileu. Atualíssima, veja a seguir:
12. A EuroSciPy 2016 será na Alemanha e a de 2015 foi na Inglaterra.
13. A SciPy Latin America 2016 será em Florianópolis agora de 16 à 20 de maio. Não perca essa chance!
14. Este é o “site” da SciPy.
15. Clicando no ícone da SciPy Library, chegamos aqui. Agora escolha Documentation, Loirinha!
16. Surfista, em Documentation escolha o Guia de Referência da SciPy.
17. Todos os tópicos tratados no Tutorial, na última versão da SciPy:
18. Em Introduction temos uma breve descrição da SciPy. Lembrem-se, os Tutoriais são para novatos!
19. É importante assinalar as convenções de importação utilizadas por padrão.
20. E agora todos os tópicos tratados na Referência, na última versão da SciPy: Assinalei os sub- pacotes que estudaremos.
21. Embarquem logo, Surfista, Loirinha e Cabelos de fogo. Vamos direto para Hogwarts.
22. Precisaremos pesquisar muito para descobrir tudo sobre Álgebra Linear Computacional na biblioteca do castelo.
23. Vou buscar com Google, Galileu! Surfista, antes disso, para chegarmos preparados, faça uma busca na Internet por Álgebra Linear Computacional.
24. Quando busquei por Numerical Linear Algebra, o 1º resultado foi da Wikipedia:
25. Não deixem de olhar as listas indicadas de software de análise numérica e de bibliotecas numéricas. Muita indicações importantes nesta dica!
26. As primeiras estações da ferrovia, antes de Hogwarts, são paradas obrigatórias, chamadas BLAS, LAPACK e NETLIB. Temos que desembarcar para conhecê-las.
27. A Estações da Luz e a Júlio Prestes, na cidade de São Paulo. Hoje também são a sede do Museu da Língua Portuguesa e da Orquestra Sinfônica do Estado de São Paulo, OSESP.
29. As bibliotecas fundamentais sobre Álgebra Linear Computacional são a BLAS e LAPACK.
30. • BLAS: Basic Linear Algebra Subprograms, • LAPACK: Linear Algebra Package. Repetindo: BLAS e LAPACK, bibliotecas altamente otimizadas que implementam os algoritmos básicos de Álgebra Linear Computacional.
31. Fui buscar BLAS no Google! Assinalei o endereço da netlib!
32. Este é o principal repositório de software científico que a comunidade científica mundial dispõe:
33. Respondendo antes de você perguntar, Loirinha:
34. A BLAS na www.netlib.org
35. FAQ é comigo! Inclusive já marquei minha 1ª pergunta:
36. Eis sua resposta, Loirinha:
37. O endereço para chegar aqui: http://www.netlib.org/blas/ No final, os links para todas as rotinas da BLAS!
38. Sim, Surfista. Clicando num link, somos encaminhados para a rotina correspondente na LAPACK. O código está em FORTRAN.
39. O código, criado pelo poderosíssimo feiticeiro Jack Dongarra em 1978. A última atualização é de 2011.
40. A LAPACK na Wikipedia:
41. A continuação da LAPACK na Wikipedia:
42. O site oficial da LAPACK, claro: www.netlib.org.
43. Agora vamos à biblioteca do Castelo de Hogwarts. É lá que encontraremos o Tutorial da scipy.linalg.
44. Marquei o motivo para usar a linalg da SciPy. Mas já tem linalg na NumPy. Por quê repetir?
45. Antes que você pergunte, Surfista, vá ler o resto desta discussão lá no Tutorial!
46. A estrutura do tutorial Linear Algebra é a seguinte: Álgebra linear Apresentação Rotinas básicas Decomposições Funções matriciais Matrizes especiais
47. Os tópicos em rotinas básicas são : Rotinas básicas Achando a inversa Resolvendo sistemas lineares Achando o determinante Calculando normas Resolvendo problemas lineares de mínimos quadrados e pseudo-inversas Inversa generalizada Vamos começar pelos determinantes.
48. O texto é um resumo da “Leitura adicional” ao final deste con junto de transparências.
49. Indo ao Reference Guide do pacote Linear Algebra, em Basics, vamos encontrar os detalhes para usar a função. Os Mestres já provaram que o cálculo de determinantes é inviável, na prática. Então, como o pacote scipy faz isto?
50. E a resposta à tua pergunta está ao pé da página Loirinha: por outro algoritmo. A utilização da função det( ) é trivial.
51. Os Mestres já ensinaram a fatoração PLU na aula passada! Insisto: Ela algo tão importante e fundamental em Álgebra linear computacional como respirar é para a vida!
52. Neste programa entramos com uma matriz A, via teclado, e em seguida chamamos a função det(A).
53. Voltando às rotinas básicas da Scipy, vejamos o que há no tópico “Calculando normas”. Rotinas básicas Achando a inversa Resolvendo sistemas lineares Achando o determinante Calculando normas Resolvendo problemas lineares de mínimos quadrados e pseudo-inversas Inversa generalizada
54. A SciPy possibilita calcular diversas normas. Tanto normas de vetores como de matrizes.
55. Por enquanto, vamos examinar apenas as normas vetoriais. Veremos as normas matriciais mais adiante. Para elas precisaremos de um pouco mais de teoria. É importantíssimo observar que todas aquelas cujo parâmetro ord é negativo não são normas.
56. Surfista, apresente alguns contra-exemplos que justifiquem a última afirmação do Mestre! Quando ord = inf temos a norma do máximo. Para ord = 1 temos a norma da soma e para ord = 2 temos norma euclidiana. Elas são anotadas ∞ , 1 e 2 respectivamente. As outras são anotadas , ( = ).
57. Já provamos que 2 = 1 2 + 2 2 + ⋯ + 2 define uma norma. A verificação que 1 = 1 + 2 + ⋯ + é uma norma é óbvia. Da mesma forma, a comprovar que ∞ = max{ 1 , 2 , ⋯ , } define uma norma também é fácil.
58. Para prová-la precisamos da desigualdade de Hölder: , ≤ ′ com ′ = −1 . A dificuldade para provar que = |1 | + |2 | + ⋯ + | | é uma norma está na desigualdade triangular: + ≤ + . Ela é conhecida na literatura como desigualdade de Minkowski.
59. Observem que a desigualdade de Hölder desempenha um papel análogo ao da desigualdade de Cauchy-Schwarz. Não vamos provar nenhuma das duas. Porém, se você estiver interessado, Surfista, já sabe as palavras-chave para buscar as demonstrações.
60. Na referência do SciPy encontramos os detalhes dos parâmetros para usar a função norma:
61. Por enquanto só veremos normas de vetores. Atenção para a referência!
62. Um exemplo de utilização da função norm( ) para vetores.
63. Mestre, parece que os valores numéricos das normas estão tendendo ao valor da ∞ ! É verdade minha filha. Prove que para ∀ ∈ ℝ temos lim →∞ = ∞
64. 1 y x A bola unitária no espaço euclidiano ℝ é o conjunto 2 = ∈ ℝ | 2 ≤ 1 . Ela descreve o conjunto 2 + 2 + 2 ≤ 1, que nossa intuição entende por uma bola (de futebol) no espaço euclidiano ℝ3. No plano seria como um CD (dos Beatles).
65. 1 1 -1 -1 y x A bola unitária no espaço ℝ com a norma da soma é o conjunto 1 = ∈ ℝ | 1 ≤ 1 . No ℝ2 ela corresponde ao conjunto definido pela desigualdade + ≤ 1, mostrado na figura.
66. 1 1 -1 -1 y x A bola unitária no espaço ℝ com a norma do máximo é o conjunto ∞ = ∈ ℝ | ∞ ≤ 1 . No ℝ2 ela corresponde ao conjunto definido pela desigualdade max{ , } ≤ 1, mostrado na figura.
67. Surfista, não vá jogar bola com 1 ou ∞ . Você pode se machucar! É, são bolas que escapam da nossa intuição euclidiana.
68. Este livro é o clássico dos clássicos sobre álgebra linear computacional. Vamos até reproduzir seu Sumário.
70. Este é outro clássico peso-pesado em álgebra linear computacional. Também vamos reproduzir seu Sumário.
71. São livros de cabeceira para quem for fazer pós-graduação em computação científica, matemática aplicada ou análise numérica.
72. O conteúdo desses dois livros é a definição mais apropriada de Álgebra linear computacional. Invista em você. Adquira os dois livros!
73. Um pouco mais de teoria. Vamos recordar o conceito de transformação linear.
74. O conceito de função é, com certeza, o mais importante da matemática. Elas são o objeto do Cálculo diferencial e integral. Em particular, vamos considerar funções de um espaço vetorial U em um outro espaço vetorial V : : → , ∈ ⟼ = () ∈ .
75. Em outras palavras, funções : → tais que: • + ⟼ ( + ) = + , • ⟼ = (). Tais funções recebem o nome especial de transformações lineares. Particularizando ainda mais: funções : → de um espaço vetorial U para outro V que preservam as operações nativas de U e V.
76. = A ℝ ℝ Matrizes × definem transformações lineares de ℝ em ℝ através da multiplicação matriz x vetor. Uma matriz 3 x 4 define uma transformação linear ∶ ℝ4 → ℝ3, através da multiplicação: 1 2 3 = 11 12 13 14 21 22 23 24 31 32 33 34 1 2 3 4
77. Lembrem-se que: Se A é uma matriz m x n então: • + = + , ∀, ∈ ℝ • = , ∀ ∈ ℝ , ∀ ∈ ℝ Que são as condições de linearidade.
78. Um exemplo simples de ℝ2 em ℝ2 é 1 2 ⟼ 1 2 = −1 2 2 −3 1 2 que corresponde ao par de igualdades ቊ 1 = −1 + 22 2 = 21 − 32
79. As equações correspondentes são ቊ 1 = 1 2 = −2 . Vemos que ela mantém a 1ª coordenada e troca o sinal da 2ª coordenada. Em outras palavras, realiza uma reflexão no eixo-x. A matriz = 1 0 0 −1 aplicada num vetor X fornece : 1 2 ⟼ 1 2 = 1 0 0 −1 1 2
80. O programa a seguir mostra a ação de uma matriz M sobre um vetor . Ele permite escolher a matriz M e o vetor .
81. Este é o programa.
82. Uma execução do programa. Em azul temos o vetor e em vermelho o vetor = X Vemos claramente que M realiza uma reflexão no eixo-x.
83. Outra execução do programa. Em azul o vetor e em vermelho o vetor = .
84. O programa da próxima transparência permitirá a visualização da ação de uma matriz M qualquer sobre um triângulo ABC escolhido livremente. Mostraremos alguns exemplos com diferentes matrizes M. Surfista, faça outros exemplos exploratórios.
85. O programa a seguir mostra a ação de uma matriz M de ordem 2 sobre um triângulo ABC.
86. Comecemos conferindo que a matriz = 1 0 0 −1 do exemplo anterior efetua uma reflexão entorno do eixo-x.
87. Entretanto a vizualização da ação da matriz = 0.5 0 0 1.5 sobre o triângulo azul não conduz a conclusões evidentes.
88. Já usando um triângulo com um lado paralelo ao eixo-x e outro paralelo ao eixo- y percebemos que M efetua uma contração no eixo-x e uma dilatação no eixo-y.
89. Já a ação desta matriz M é complicada de descrever.
90. A matriz M definida por = cos() −() () cos() , quando aplicada num vetor v efetua uma rotação de ângulo no sentido horário em v. Isto é, se = então o ângulo entre v e u, medido de v para u é de radianos.
91. Fiz um programa que recebe um quadrado arbitrariamente escolhido e gira-o no sentido horário de uma quantidade escolhida de graus..
92. Uma rotação de 45º e outra de 90º. Notem que na rotação a norma dos vetores (e portanto dos objetos construídos com eles) permanece inalterada.
93. O código do programa. Ficou longo porque optei pela didática.
94. O final do código do programa..
95. x x I A matriz identidade I é uma matriz quadrada de ordem n com 1’s na diagonal e 0’s fora dela. A identidade, como transformação linear não faz nada! ℝ ℝ I x = x
96. Uma matriz M muito próxima da identidade = 1 0 0 1 . Se M fosse a identidade I veríamos apenas o triângulo rosa, cobrindo o azul.
97. x y A z B ℝ ℝ ℝ BA O produto BA de uma matriz B por uma matriz A corresponde à aplicação linear composta de A seguida por B:
98. Lembrem-se que a multiplicação de matrizes não é comutativa. Na aula passada exibimos um programa que mostra que, usualmente, ≠ . Por exemplo, para = 1 2 1 2 e = −1 1 1 −1 temos = 1 −1 1 −1 e = 0 0 0 0 .
99. A inversa de uma matriz quadrada A de ordem n é uma matriz quadrada B de ordem n tal que ∙ = ∙ = Existem matrizes que não possuem inversa – são não-inversíveis. Uma matriz não possui mais que uma inversa. Além disso, a inversa da inversa é a própria!
100. Prova-se que A é inversível quando, e apenas quando, é não-singular, i.é, det() ≠ 0 Entretanto, verificar se det() ≠ 0 para saber se A é inversível é uma técnica nunca utilizada em Álgebra linear computacional.
101. x y A A-1 ℝ ℝ x y A x A-1 ℝ ℝ ℝ I Pela própria definição, a inversa A-1 de uma matriz A , quando considerada como aplicação linear, desfaz a ação de A.
102. Voltando às rotinas básicas do foquemos nossa atenção no tópico “Achando a inversa”. Rotinas básicas Achando a inversa Resolvendo sistemas lineares Achando o determinante Calculando normas Resolvendo problemas lineares de mínimos quadrados e pseudo-inversas Inversa generalizada
103. Vejam o que Tutorial do scipy.linalg informa sobre a inversa de uma matriz:
104. Indo à Referência e clicando em inv( ) obtive as informações abaixo:
105. Obtendo a inversa de uma matriz A e conferindo.
106. Vejam a ação de uma matriz não-inversível M sobre um triângulo: Não há como desfazer a ação da M sobre o trtiângulo!
107. Essa matriz é quase não- inversível. Vejam como ela achata um triângulo:
108. Pedi ao Mestre para fazer um programa parecido àquele que mostra o efeito visual da ação de uma matriz num triângulo. Mas para a matriz inversa! Ele fará Loirinha, você é a queridinha dele!
109. O seu pedido satisfeito, Loirinha!
110. Verificando as ações de uma matriz e sua inversa, como transformações lineares. O deslocamento é um artifício para enxergarmos o efeito da inversa (em verde).
111. Toda esta recordação sobre o efeito sobre polígonos das transformações lineares associadas às matrizes será utilizada na próxima aula!
112. Tchau, até a próxima aula!