#!/usr/bin/perl -w
use Math::GSL          qw/:all/;
use Math::GSL::Eigen   qw/:all/;
use Math::GSL::Matrix  qw/:all/;
use Math::GSL::Vector  qw/:all/;
use Math::GSL::Complex qw/:all/;
use Math::GSL::Errno qw/:all/;
use Math::Complex;
use Data::Dumper;
use strict;
my $matrix = Math::GSL::Matrix->new(2,2)
                              ->set_row(0, [0,-1] )
                              ->set_row(1, [1, 0] );
print <<STUFF;
Finding eigenvalue/eigenvectors for
[ 0  -1 ]
[ 1   0 ]
STUFF

my $evec   = gsl_matrix_complex_alloc(2,2);
my $eigen  = gsl_eigen_nonsymmv_alloc(2);
my $vector = gsl_vector_complex_alloc(2);

# this actually computes the eigenvalues and vectors
gsl_eigen_nonsymmv($matrix->raw,$vector, $evec, $eigen);

my $x = gsl_vector_complex_real($vector);
my $y = gsl_vector_complex_imag($vector);

my $eig1 = as_complex( gsl_vector_get($x->{vector}, 0) ,  gsl_vector_get($y->{vector}, 0) );
my $eig2 = as_complex( gsl_vector_get($x->{vector}, 1) ,  gsl_vector_get($y->{vector}, 1) );

my $u1 = as_complex( gsl_matrix_complex_get($evec, 0, 0) );
my $u2 = as_complex( gsl_matrix_complex_get($evec, 1, 0) );
my $v1 = as_complex( gsl_matrix_complex_get($evec, 0, 1) );
my $v2 = as_complex( gsl_matrix_complex_get($evec, 1, 1) );

print <<ANSWER;

Eigenvectors:

u = ($u1,$u2)
v = ($v1,$v2)

Eigenvalues:

lambda_0 = $eig1
lambda_1 = $eig2

ANSWER

sub as_complex {
    my ($w,$v) = @_;
    my ($x,$y);
    if( ref $w ) {
        ($x,$y) = (gsl_real($w),gsl_imag($w));
    } else {
        ($x,$y) = ($w,$v);
    }
    my $z      = Math::Complex->make( $x, $y);
    return qq{$z};
}