15 Autosistemas

Esse capítulo descreve funções para cálculo de autovalores e autovetores de matrizes. Existem rotinas para autosistemas de diversos tipos: simétrico real, não simétrico real, de Hermite complexo, definido simétrico generalizado real, de Hermite definido generalizado complexo, e não simétrico generalizado real. Autovalores podem ser calculados com ou sem autovetores. Os algoritmos de matriz de Hermite real e simétrica são bidiagonalização simétrica seguido por redução QR. O algoritmo não simétrico é o QR de Francis dupla jornada. O algoritmo não simétrico generalizado é o método QZ devido a Moler e Stewart.

As funções descritas nesse capítulo são declaradas no arquivo de cabeçalho ‘gsl_eigen.h’.

15.1 Matrizes Simétricas Reais

Para matrizes simétricas e reais, a biblioteca usa a bidiagonalização simétrica e o método de redução QR. Isso é descrito em Golub & van Loan, seção 8.3. Os autovalores calculados são precisos para uma exatidão absoluta de \epsilon ||A||_2, onde \epsilon é a precisão da máquina.

Function: gsl_eigen_symm_workspace * gsl_eigen_symm_alloc (const size_t n)

Essa função aloca um espaço de trabalho para o cálculo de autovalores das matrizes simétricas e reais de ordem n-por-n. O tamanho do espaço de trabalho é O(2n) (39).

Function: void gsl_eigen_symm_free (gsl_eigen_symm_workspace * w)

Essa função libera a memória associada ao espaço de trabalho w.

Function: int gsl_eigen_symm (gsl_matrix * A, gsl_vector * eval, gsl_eigen_symm_workspace * w)

Essa função calcula os autovalores da matriz simétrica real A. Espaço de trabalho adicional de tamanho apropriado deve ser fornecido em w. A diagonal e a parte triangular baixa de A são destruidas durante o cálculo, mas a parte triangular estritamente alta não é referenciada. Os autovalores são armazenados no vetor eval e são desordenados.

Function: gsl_eigen_symmv_workspace * gsl_eigen_symmv_alloc (const size_t n)

Essa função aloca um espaço de trabalho para o cálculo de autovalores e de autovetores de matrizes simétricas e reais de ordem n-por-n . O tamanho do espaço de trabalho é O(4n).

Function: void gsl_eigen_symmv_free (gsl_eigen_symmv_workspace * w)

Essa função libera a memória associada ao espaço de trabalho w.

Function: int gsl_eigen_symmv (gsl_matrix * A, gsl_vector * eval, gsl_matrix * evec, gsl_eigen_symmv_workspace * w)

Essa função calcula os autovalores e os autovetores da matriz simétrica real A. Espaço de trabalho adicional de tamanho apropriado deve ser fornecido em w. A diagonal e a parte triangular baixa de A são destruidas durante a computação, mas a parte triangular estritamente alta não é referenciada. Os autovalores são armazenados no vetor eval e são desordenados. Os correspondentes autovetores são armazenados em colunas da matrix evec. Por exemplo, o autovetor na primeira coluna corresponde ao primeiro autovalor. Os autovetores são garantidamente serem mutuamente ortogonais e normalizados para módulo unitário.

15.2 Matrizes de Hermite Complexas

Para matrizes de Hermite, a biblioteca usa a forma complexa da bidiagonalização simétrica e o método de redução QR.

Function: gsl_eigen_herm_workspace * gsl_eigen_herm_alloc (const size_t n)

Essa função aloca um espaço de trabalho para o cálculo de autovalores de matrizes de Hermite complexas de ordem n-por-n. O tamanho do espaço de trabalho é O(3n).

Function: void gsl_eigen_herm_free (gsl_eigen_herm_workspace * w)

Essa função libera a memória associada com o espaço de trabalho w.

Function: int gsl_eigen_herm (gsl_matrix_complex * A, gsl_vector * eval, gsl_eigen_herm_workspace * w)

Essa função calcula os autovalores da matriz de Hermite complexa A. Espaço de trabalho adicional de tamanho apropriado deve ser fornecido em w. A diagonal e a parte triangular baixa de A são destruídos durante o cálculo, mas a parte triangular estritamente alta não é referenciada. As partes imaginárias da diagonal são assumidas serem zero e não são referenciadas. Os autovalores são armazenados no vetor eval e são desordenadas.

Function: gsl_eigen_hermv_workspace * gsl_eigen_hermv_alloc (const size_t n)

Essa função aloca um espaço de trabalho para o cálculo de autovalores e autovetores de matrizes de Hermite complexas de ordem n-por-n. O tamanho do espaço de trabalho é O(5n).

Function: void gsl_eigen_hermv_free (gsl_eigen_hermv_workspace * w)

Essa função libera a memória associada ao espaço de trabalho w.

Function: int gsl_eigen_hermv (gsl_matrix_complex * A, gsl_vector * eval, gsl_matrix_complex * evec, gsl_eigen_hermv_workspace * w)

Essa função calcula os autovalores e autovetores da matriz de Hermite complexa A. Espaço de trabalho adicional de tamanho apropriado deve ser fornecido em w. A diagonal e a parte triangular baixa de A são destruídos durante o cálculo, mas a parte triangular estritamente alta não é referenciada. As partes imaginárias da diagonal são assumidas serem zero e não são referenciadas. Os autovalores são armazenados no vetor eval e são desordenados. Os correspondentes autovetores complexos são armazenados nas colunas da matrix evec. Por exemplo, o autovetor na primeira coluna corresponde ao primeiro autovalor. Os autovetores são garantidamente serem mutualmente ortogonal e normalizados para módulo unitário.

15.3 Matrizes Assimétricas Reais

A solução do problema de autosistema não simétrico real para uma matriz A envolve cálculos da decomposição de Schur

Function: gsl_eigen_nonsymm_workspace * gsl_eigen_nonsymm_alloc (const size_t n)

Essa função aloca um espaço de trabalho para o cálculo de autovalores de matrizes não simétricas reais de ordem n-por-n. O tamanho do espaço de trabalho é O(2n).

Function: void gsl_eigen_nonsymm_free (gsl_eigen_nonsymm_workspace * w)

Essa função libera a memória associada ao espaço de trabalho w.

Function: void gsl_eigen_nonsymm_params (const int compute_t, const int balance, gsl_eigen_nonsymm_workspace * w)

Essa função ajusta alguns parâmetros que determinam como o problema de autovalor é resolvido em chamadas subsequentes a gsl_eigen_nonsymm.

Se compute_t é ajustada para 1, a forma completa de Schur T irá ser calculada por gsl_eigen_nonsymm. Se compute_t for ajustada para 0, T, a forma completa de Schur, não será calculada (esse é o ajuste padrão). Calculando a forma completa de Schur T requer aproximadamente 1.5–2 vezes o número de flops (40).

Se balance for ajustada para 1, uma transformação de balanceamento é aplicada à matriz principal para o cálculo de autovalores. Essa transoformação tem por objetivo fazer com que as linhas e colunas da matriz tenham normas comparáveis, e pode resultar em autovalores mais precisos para matrizes cujas entradas variam grandemente em módulo. Veja Balanceamento para mais informação. Note que a transformação de balanceamento não preserva a ortogonalidade dos vetores de Schur, de forma que se você deseja calcular os vetores de Schur com gsl_eigen_nonsymm_Z você irá obter os vetores de Schur da matriz balanceada ao invés da matriz original. A relação irá ser

T = Qt D−1 A D Q

onde Q é a matriz dos vetores de Schur para a matriz balanceada, e D é a transformação de balanceamento. Então gsl_eigen_nonsymm_Z irá calcular uma matriz Z que satisfaz

T = Z−1 A Z

com Z = D Q. Note que Z não irá ser orthogonal. Por essa razão, o balanceamento não é executado por padrão.

Function: int gsl_eigen_nonsymm (gsl_matrix * A, gsl_vector_complex * eval, gsl_eigen_nonsymm_workspace * w)

Essa função calcula os autovalores da matriz não simétrica real A e armazena os autovalores no vetor eval. Se T for desejada, é armazenada na porção alta de A na saída. De outra forma, na saída, a diagonal de A irá conter os autovalores reais de tamanho 1-por-1 e sistemas de autovalores conjugados complexos de tamanho 2-por-2, e o restante de A é destruído. Em raros casos, essa função pode falhar em encontrar todos os autovalores. Se isso ocorrer, um código de erro é retornado e o número de autovalores convergidos é armazenado em w->n_evals. Os autovalores convergidos são armazenados no início de eval.

Function: int gsl_eigen_nonsymm_Z (gsl_matrix * A, gsl_vector_complex * eval, gsl_matrix * Z, gsl_eigen_nonsymm_workspace * w)

Essa função é identica a gsl_eigen_nonsymm exceto que gsl_eigen_nonsymm_Z também calcula os vetores de Schur e armazena-os em Z.

Function: gsl_eigen_nonsymmv_workspace * gsl_eigen_nonsymmv_alloc (const size_t n)

Essa função aloca um espaço de trabalho para o cálculo de autovalores e autovetores de matrizes não simétricas reais de ordem n-por-n. O tamanho do espaço de trabalho é O(5n).

Function: void gsl_eigen_nonsymmv_free (gsl_eigen_nonsymmv_workspace * w)

Essa função libera a memória associada ao espaço de trabalho w.

Function: void gsl_eigen_nonsymmv_params (const int balance, gsl_eigen_nonsymm_workspace * w)

Essa função ajusta parâmetros que determinam como o problema de autovalor é resolvido em subsequêntes chamadas a gsl_eigen_nonsymmv. Se balance for ajustado para 1, uma transformação de balanceamento é aplicada à matriz. Veja gsl_eigen_nonsymm_params para mais informação. O Balanceamento é desabilitado por padrão uma vez que não preserva a ortogonalidade nos vetores de Schur.

Function: int gsl_eigen_nonsymmv (gsl_matrix * A, gsl_vector_complex * eval, gsl_matrix_complex * evec, gsl_eigen_nonsymmv_workspace * w)

Essa função calcula autovalores e autovetores direitos da matriz não simétrica real A de ordem n-por-n. Essa função primeiramente chama gsl_eigen_nonsymm para calcular os autovalores, a forma de Schur T, e vetores de Schur. Então encontra dos autovetores de T e transforma-os de volta usando os vetores de Schur. Os vetores de Schur são destruídos no processo, mas podem ser gravados usando gsl_eigen_nonsymmv_Z. Os autovetores calculados são normalizados para terem módulo unitário. Na saída, a porção alta de A contém a forma de Schur T. Se gsl_eigen_nonsymm vier a falhar, nenhum autovetor é calculado, e um código de erro é retornado.

Function: int gsl_eigen_nonsymmv_Z (gsl_matrix * A, gsl_vector_complex * eval, gsl_matrix_complex * evec, gsl_matrix * Z, gsl_eigen_nonsymmv_workspace * w)

Essa função é idêntica a gsl_eigen_nonsymmv exceto que gsl_eigen_nonsymmv_Z também grava os vetores de Schur em Z.

15.4 Autosistemas Definidos Simétricos Generalizados Reais

O problema de autovalor definido simétrico generalizado real é encontrar autovalores \lambda e autovetores x tais que

Function: gsl_eigen_gensymm_workspace * gsl_eigen_gensymm_alloc (const size_t n)

Essa função aloca um espaço de trabalho para o cálculo de autovalores de autosistemas definidos simétricos generalizados reais de ordem n-por-n. O tamanho do espaço de trabalho é O(2n).

Function: void gsl_eigen_gensymm_free (gsl_eigen_gensymm_workspace * w)

Essa função libera a memória associada ao espaço de trabalho w.

Function: int gsl_eigen_gensymm (gsl_matrix * A, gsl_matrix * B, gsl_vector * eval, gsl_eigen_gensymm_workspace * w)

Essa função calcula os autovalores do par de matrizes definidas simétricas generalizadas reais (A, B), e armazena-os em eval, usando o método mostrado acima. Na saída, B contém sua decomposição de Cholesky e A é destruída.

Function: gsl_eigen_gensymmv_workspace * gsl_eigen_gensymmv_alloc (const size_t n)

Essa função aloca um espaço de trabalho para o cálculo de autovalores e autovetores de autosistemas definidos simétricos generalizados reais de ordem n-por-n. O tamanho do espaço de trabalho é O(4n).

Function: void gsl_eigen_gensymmv_free (gsl_eigen_gensymmv_workspace * w)

Essa função libera a memária associada ao espaço de trabalho w.

Function: int gsl_eigen_gensymmv (gsl_matrix * A, gsl_matrix * B, gsl_vector * eval, gsl_matrix * evec, gsl_eigen_gensymmv_workspace * w)

Essa função calcula os autovalores e autovetores do par de matrizes definidas simétricas generalizadas reais (A, B), e armazena-as em eval e em evec respectivamente. Os autovetores calculados são normalizados para terem módulo unitário. Na saída, B contém sua decomposição de Cholesky e A é destruída.

15.5 Autosistemas de Hermite Definidos Generalizados Complexos

O problema do autovalor definido de Hermite generalizado complexo é encontrar autovalores \lambda e autovetores x tais que

Function: gsl_eigen_genherm_workspace * gsl_eigen_genherm_alloc (const size_t n)

Essa função aloca um espaço de trabalho para o cálculo de autovalores do autosistema definido de Hermite generalizado complexo de ordem n-por-n. O tamanho do espaço de trabalho é O(3n).

Function: void gsl_eigen_genherm_free (gsl_eigen_genherm_workspace * w)

Essa função libera a memória associada ao espaço de trabalho w.

Function: int gsl_eigen_genherm (gsl_matrix_complex * A, gsl_matrix_complex * B, gsl_vector * eval, gsl_eigen_genherm_workspace * w)

Essa função calcula os autovalores do par de matrizes definidas de Hermite generalizadas complexas (A, B), e armazena-os em eval, usando o método mostrado acima. Na saída, B contém sua decomposição de Cholesky e A é destruída.

Function: gsl_eigen_genhermv_workspace * gsl_eigen_genhermv_alloc (const size_t n)

Essa função aloca um espaço de trabalho para o cálculo de autovalores e de autovetores de autosistemas definidos de Hermite generalizados complexos de ordem n-por-n. O tamanho do espaço de trabalho é O(5n).

Function: void gsl_eigen_genhermv_free (gsl_eigen_genhermv_workspace * w)

Essa função libera a memória associada ao espaço de trabalho w.

Function: int gsl_eigen_genhermv (gsl_matrix_complex * A, gsl_matrix_complex * B, gsl_vector * eval, gsl_matrix_complex * evec, gsl_eigen_genhermv_workspace * w)

Essa função calcula os autovalores e autovetores do par de matrizes definidas de Hermite generalizadas complexas (A, B), e armazena-os em eval e evec respectivamente. Os autovetores calculados são normalizados para terem módulo unitário. Na saída, B contém sua decomposição de Cholesky e A é destruída.

15.6 Autosistemas Não Simétricos Generalizados Reais

Fornecidas duas matrizes quadradas (A, B), o problema de autovalor não simétrico generalizado é encontrar autovalores \lambda e autovetores x tais que

Se o determinante da matriz pencil A - \lambda B for zero para todos os \lambda, o problema é dito singular; de outra forma o problema é chamado regular. Singularidade normalmente leva a algum \alpha = \beta = 0 o que significa que o autoproblema é de condicionamento hostil e geralmente não tem soluções de autovalores bem definidas. As rotinas abaixo foram pensadas para matrizes pencils regulares e podem retornar resultados imprevisíveis quando aplicadas a pencils singulares.

A solução do problema de autosistema não simétrico generalizado real para um par de matrizes (A, B) envolve cálculos da decomposição de Schur generalizada

Function: gsl_eigen_gen_workspace * gsl_eigen_gen_alloc (const size_t n)

Essa função aloca um espaço de trabalho para o cálculo de autovalores de autosistemas não simétricos generalizados reais de ordem n-por-n. O tamanho do espaço de trabalho é O(n).

Function: void gsl_eigen_gen_free (gsl_eigen_gen_workspace * w)

Essa função libera a memória associada ao espaço de trabalho w.

Function: void gsl_eigen_gen_params (const int compute_s, const int compute_t, const int balance, gsl_eigen_gen_workspace * w)

Essa função ajusta alguns parâmetros que determinam como o problema de autovalor é resolvido em subsequêntes chamadas a gsl_eigen_gen.

Se compute_s for ajustado para 1, a forma completa de Schur S irá ser calculada por gsl_eigen_gen. Se compute_s for ajustada para 0, S não irá ser calculada (esse é o ajuste padrão). S é uma matriz quase triangular alta com blocos de tamanho 1-por-1 e 2-por-2 sobre sua diagonal. Blocos 1-por-1 correspondem a autovalores reais, e blocos 2-por-2 correspondem a autovalores complexos.

Se compute_t for ajustada para 1, a forma completa de Schur T irá ser calculada por gsl_eigen_gen. Se compute_t for ajustada para 0, T não irá ser calculada (esse é o ajuste padrão). T é uma matriz triangular alta com elementos não negativos em sua diagonal. Quaisquer blocos 2-por-2 em S irá corresponder a um bloco diagonal 2-por-2 em T.

O parâmetro balance é atualmetne ignorado, uma vez que balanceamento generalizado não está ainda implementado.

Function: int gsl_eigen_gen (gsl_matrix * A, gsl_matrix * B, gsl_vector_complex * alpha, gsl_vector * beta, gsl_eigen_gen_workspace * w)

Essa função calcula os autovalores do par de matrizes não simétrico generalizado real (A, B), e armazena-os como pares em (alpha, beta), onde alpha é complexo e beta é real. Se \beta_i for não nulo, então \lambda = \alpha_i / \beta_i é um autovalor. Da mesma forma, se \alpha_i for não nulo, então \mu = \beta_i / \alpha_i é um autovalor do problema alternativo \mu A y = B y. Os elementos de beta são normalizados para serem não negativos.

Se S for desejado, é armazenado em A na saída. Se T for desejado, é armazenado em B na saída. A ordem de autovalores em (alpha, beta) segue a ordem dos blocos de diagonal nas formas de Schur S e T. Em raros casos, essa função pode falhar em encontrar todos os autovalores. Se isso ocorrer, um código de erro é retornado.

Function: int gsl_eigen_gen_QZ (gsl_matrix * A, gsl_matrix * B, gsl_vector_complex * alpha, gsl_vector * beta, gsl_matrix * Q, gsl_matrix * Z, gsl_eigen_gen_workspace * w)

Essa função é idêntica a gsl_eigen_gen exceto que também calcula os vetores de Schur esquerdo e direito e armazena-os em Q e Z respectivamente.

Function: gsl_eigen_genv_workspace * gsl_eigen_genv_alloc (const size_t n)

Essa função aloca um espaço de trabalho para o cálculo de autovalores e autovetores de autosistemas não simétricos generalizados de ordem n-por-n. O tamanho do espaço de trabalho é O(7n).

Function: void gsl_eigen_genv_free (gsl_eigen_genv_workspace * w)

Essa função libera a memória associada ao espaço de trabalho w.

Function: int gsl_eigen_genv (gsl_matrix * A, gsl_matrix * B, gsl_vector_complex * alpha, gsl_vector * beta, gsl_matrix_complex * evec, gsl_eigen_genv_workspace * w)

Essa função calcula autovalores e autovetores direitos do par de matrizes não simétricas generalizadas reais de ordem n-por-n (A, B). Os autovalores são armazenados em (alpha, beta) e os autovetores são armazenados em evec. Essa função primeiro chama gsl_eigen_gen para calcular os autovalores, formas de Schur, e vetores de Schur. Então encontra autovetores das formas de Schur e transforma-as de volta usando os vetores de Schur. Os vetores de Schur são destruídos no processo, mas podem ser salvos usando gsl_eigen_genv_QZ. Os autovetores calculados são normalizados para terem módulo unitário. Na saída, (A, B) possuem a forma generalizada de Schur (S, T). Se gsl_eigen_gen falha, nenhum autovetor é calculado, e um código de erro é retornado.

Function: int gsl_eigen_genv_QZ (gsl_matrix * A, gsl_matrix * B, gsl_vector_complex * alpha, gsl_vector * beta, gsl_matrix_complex * evec, gsl_matrix * Q, gsl_matrix * Z, gsl_eigen_genv_workspace * w)

Essa função é idêntica a gsl_eigen_genv exceto que também calcula os vetores de Schur esquerdo e direito e armazena-os em Q e Z respectivamente.

15.7 Ordenando Autovalores e Autovetores

Function: int gsl_eigen_symmv_sort (gsl_vector * eval, gsl_matrix * evec, gsl_eigen_sort_t sort_type)

Essa função simultâneamente ordena os autovalores armazenados no vetor eval e os correspondentes autovetores reais armazenados nas colunas da matriz evec em ordem ascendente ou descendente conforme o valor do parâmetro sort_type,

GSL_EIGEN_SORT_VAL_ASC

ordem ascendente em valor numérico

GSL_EIGEN_SORT_VAL_DESC

ordem descendente em valor numérico

GSL_EIGEN_SORT_ABS_ASC

ordem ascendente de módulo

GSL_EIGEN_SORT_ABS_DESC

ordem descendente de módulo

Function: int gsl_eigen_hermv_sort (gsl_vector * eval, gsl_matrix_complex * evec, gsl_eigen_sort_t sort_type)

Essa função ordena simultâneamente os autovalores armazenados no vetor eval e os correspondentes autovetores complexos armazenados nas colunas da matriz evec em ordem ascendente ou descendente conforme o valor deo parâmetro sort_type como mostrado acima.

Function: int gsl_eigen_nonsymmv_sort (gsl_vector_complex * eval, gsl_matrix_complex * evec, gsl_eigen_sort_t sort_type)

Essa função ordena simultâneamente os autovalores armazenados no vetor eval e os correspondentes autovetores complexos armazenados nas colunas da matrix evec em ordem ascendente ou descendente conforme o valor do parâmetro sort_type como mostrado acima. Somente GSL_EIGEN_SORT_ABS_ASC e GSL_EIGEN_SORT_ABS_DESC são suportadas devido aos autovalores serem complexos.

Function: int gsl_eigen_gensymmv_sort (gsl_vector * eval, gsl_matrix * evec, gsl_eigen_sort_t sort_type)

Essa função ordena simultâneamente os autovalores armazenados no vetor eval os correspondentes autovetores reais armazenados nas colunas da matriz evec em ordem ascendente ou descendente conforme o valor do parâmetro sort_type como mostrado acima.

Function: int gsl_eigen_genhermv_sort (gsl_vector * eval, gsl_matrix_complex * evec, gsl_eigen_sort_t sort_type)

Essa função ordena simultâneamente os autovalores armazenados no vetor eval e os correspondentes autovetores complexos armazenados nas colunas da matriz evec em ordem ascendente ou descendente conforme o valor do parâmetro sort_type como mostrado acima.

Function: int gsl_eigen_genv_sort (gsl_vector_complex * alpha, gsl_vector * beta, gsl_matrix_complex * evec, gsl_eigen_sort_t sort_type)

Essa função ordena simultâneamente os autovalores armazenados nos vetores (alpha, beta) e os correspondentes autovetores complexos armazenados nas colunas da matriz evec em ordem ascendente ou descendente conforme o valor do parâmetro sort_type como mostrado acima. Somente GSL_EIGEN_SORT_ABS_ASC e GSL_EIGEN_SORT_ABS_DESC são suportados devido aos autovalores serem complexos.

15.8 Exemplos

Os seguinte programas calculam os autovalores e autovetores da matriz de Hilbert de quarta ordem, H(i,j) = 1/(i + j + 1).

#include <stdio.h>
#include <gsl/gsl_math.h>
#include <gsl/gsl_eigen.h>

int
main (void)
{
  double data[] = { 1.0  , 1/2.0, 1/3.0, 1/4.0,
                    1/2.0, 1/3.0, 1/4.0, 1/5.0,
                    1/3.0, 1/4.0, 1/5.0, 1/6.0,
                    1/4.0, 1/5.0, 1/6.0, 1/7.0 };

  gsl_matrix_view m 
    = gsl_matrix_view_array (data, 4, 4);

  gsl_vector *eval = gsl_vector_alloc (4);
  gsl_matrix *evec = gsl_matrix_alloc (4, 4);

  gsl_eigen_symmv_workspace * w = 
    gsl_eigen_symmv_alloc (4);
  
  gsl_eigen_symmv (&m.matrix, eval, evec, w);

  gsl_eigen_symmv_free (w);

  gsl_eigen_symmv_sort (eval, evec, 
                        GSL_EIGEN_SORT_ABS_ASC);
  
  {
    int i;

    for (i = 0; i < 4; i++)
      {
        double eval_i 
           = gsl_vector_get (eval, i);
        gsl_vector_view evec_i 
           = gsl_matrix_column (evec, i);

        printf ("eigenvalue = %g\n", eval_i);
        printf ("eigenvector = \n");
        gsl_vector_fprintf (stdout, 
                            &evec_i.vector, "%g");
      }
  }

  gsl_vector_free (eval);
  gsl_matrix_free (evec);

  return 0;
}

Aqui está o início da saída do programa,

$ ./a.out 
eigenvalue = 9.67023e-05
eigenvector = 
-0.0291933
0.328712
-0.791411
0.514553
...

A saída acima pode ser comparada à correspondente saída do GNU OCTAVE,

octave> [v,d] = eig(hilb(4));
octave> diag(d)  
ans =

   9.6702e-05
   6.7383e-03
   1.6914e-01
   1.5002e+00

octave> v 
v =

   0.029193   0.179186  -0.582076   0.792608
  -0.328712  -0.741918   0.370502   0.451923
   0.791411   0.100228   0.509579   0.322416
  -0.514553   0.638283   0.514048   0.252161

Note que os autovetores podem diferir de uma mudança de sinal, uma vez que o sinal de um autovetor é arbitrário.

O seguinte programa ilustra o uso do autoresolvedor não simétrico, através do cálculo dos autovalores e dos autovetores da matriz de Vandermonde V(x;i,j) = x_i^{n - j} com x = (-1,-2,3,4).

#include <stdio.h>
#include <gsl/gsl_math.h>
#include <gsl/gsl_eigen.h>

int
main (void)
{
  double data[] = { -1.0, 1.0, -1.0, 1.0,
                    -8.0, 4.0, -2.0, 1.0,
                    27.0, 9.0, 3.0, 1.0,
                    64.0, 16.0, 4.0, 1.0 };

  gsl_matrix_view m 
    = gsl_matrix_view_array (data, 4, 4);

  gsl_vector_complex *eval = gsl_vector_complex_alloc (4);
  gsl_matrix_complex *evec = gsl_matrix_complex_alloc (4, 4);

  gsl_eigen_nonsymmv_workspace * w = 
    gsl_eigen_nonsymmv_alloc (4);
  
  gsl_eigen_nonsymmv (&m.matrix, eval, evec, w);

  gsl_eigen_nonsymmv_free (w);

  gsl_eigen_nonsymmv_sort (eval, evec, 
                           GSL_EIGEN_SORT_ABS_DESC);
  
  {
    int i, j;

    for (i = 0; i < 4; i++)
      {
        gsl_complex eval_i 
           = gsl_vector_complex_get (eval, i);
        gsl_vector_complex_view evec_i 
           = gsl_matrix_complex_column (evec, i);

        printf ("eigenvalue = %g + %gi\n",
                GSL_REAL(eval_i), GSL_IMAG(eval_i));
        printf ("eigenvector = \n");
        for (j = 0; j < 4; ++j)
          {
            gsl_complex z = 
              gsl_vector_complex_get(&evec_i.vector, j);
            printf("%g + %gi\n", GSL_REAL(z), GSL_IMAG(z));
          }
      }
  }

  gsl_vector_complex_free(eval);
  gsl_matrix_complex_free(evec);

  return 0;
}

Aqui está o começo da saída do programa,

$ ./a.out 
eigenvalue = -6.41391 + 0i
eigenvector = 
-0.0998822 + 0i
-0.111251 + 0i
0.292501 + 0i
0.944505 + 0i
eigenvalue = 5.54555 + 3.08545i
eigenvector = 
-0.043487 + -0.0076308i
0.0642377 + -0.142127i
-0.515253 + 0.0405118i
-0.840592 + -0.00148565i
...

A saída acima pode ser comparada à saída correspondente do GNU OCTAVE,

octave> [v,d] = eig(vander([-1 -2 3 4]));
octave> diag(d)
ans =

  -6.4139 + 0.0000i
   5.5456 + 3.0854i
   5.5456 - 3.0854i
   2.3228 + 0.0000i

octave> v
v =

 Columns 1 through 3:

  -0.09988 + 0.00000i  -0.04350 - 0.00755i  -0.04350 + 0.00755i
  -0.11125 + 0.00000i   0.06399 - 0.14224i   0.06399 + 0.14224i
   0.29250 + 0.00000i  -0.51518 + 0.04142i  -0.51518 - 0.04142i
   0.94451 + 0.00000i  -0.84059 + 0.00000i  -0.84059 - 0.00000i

 Column 4:

  -0.14493 + 0.00000i
   0.35660 + 0.00000i
   0.91937 + 0.00000i
   0.08118 + 0.00000i

Note que os autovetores correspondem aos autovalores 5.54555 + 3.08545i diferem através de uma constante multiplicativa 0.9999984 + 0.0017674i que é um fator arbitrário de fase de módulo 1.

15.9 Referências e Leituras Adicionais

Informação adicional sobre os algoritmos descritos nessa seção podem ser encontradas no seguinte livro,

• G. H. Golub, C. F. Van Loan, Matrix Computations (3rd Ed, 1996), Johns Hopkins University Press, ISBN 0-8018-5414-8.

Informação adicional sobre o algoritmo QZ para autosistemas generalizados podem ser encontrados no seguinte artigo,

• C. Moler, G. Stewart, “An Algorithm for Generalized Matrix Eigenvalue Problems”, SIAM J. Numer. Anal., Vol 10, No 2, 1973.

Rotinas de autosistemas para matrizes muito grandes podem ser encontrados na biblioteca Fortran LAPACK. A bilioteca LAPACK é descrita em,

• LAPACK Users’ Guide (Terceira Edição, 1999), Publicado pela SIAM, ISBN 0-89871-447-8.

  http://www.netlib.org/lapack

O código fonte da LAPACK pode ser encontrado no website acima juntamente com uma cópia online copy do guia de usuários.

Esse documento foi gerado em 23 de Julho de 2013 usando texi2html 5.0.