1 package Math
::GSL
::Eigen
::Test
;
2 use base
q{Test::Class};
4 use Math
::GSL qw
/:all/;
5 use Math
::GSL
::Test qw
/:all/;
6 use Math
::GSL
::Eigen qw
/:all/;
7 use Math
::GSL
::Matrix qw
/:all/;
8 use Math
::GSL
::Vector qw
/:all/;
9 use Math
::GSL
::Complex qw
/:all/;
10 use Math
::GSL
::Machine qw
/:all/;
14 sub make_fixture
: Test
(setup
) {
16 $self->{eigen
} = gsl_eigen_symm_alloc
(2);
19 sub teardown
: Test
(teardown
) {
22 sub GSL_EIGEN_SYMM_ALLOC
: Tests
{
24 isa_ok
($self->{eigen
}, 'Math::GSL::Eigen');
27 sub GSL_EIGEN_SYMMV_ALLOC
: Tests
{
28 my $eigen = gsl_eigen_symmv_alloc
(5);
29 isa_ok
($eigen, 'Math::GSL::Eigen');
32 sub GSL_EIGEN_HERM_ALLOC
: Tests
{
33 my $eigen = gsl_eigen_herm_alloc
(5);
34 isa_ok
($eigen, 'Math::GSL::Eigen');
37 sub GSL_EIGEN_HERMV_ALLOC
: Tests
{
38 my $eigen = gsl_eigen_hermv_alloc
(5);
39 isa_ok
($eigen, 'Math::GSL::Eigen');
42 sub GSL_EIGEN_NONSYMM_ALLOC
: Tests
{
43 my $eigen = gsl_eigen_nonsymm_alloc
(5);
44 isa_ok
($eigen, 'Math::GSL::Eigen');
47 sub GSL_EIGEN_NONSYMMV_ALLOC
: Tests
{
48 my $eigen = gsl_eigen_nonsymmv_alloc
(5);
49 isa_ok
($eigen, 'Math::GSL::Eigen');
52 sub GSL_EIGEN_GENSYMM_ALLOC
: Tests
{
53 my $eigen = gsl_eigen_gensymm_alloc
(5);
54 isa_ok
($eigen, 'Math::GSL::Eigen');
57 sub GSL_EIGEN_GENSYMMV_ALLOC
: Tests
{
58 my $eigen = gsl_eigen_gensymmv_alloc
(5);
59 isa_ok
($eigen, 'Math::GSL::Eigen');
62 sub GSL_EIGEN_GENHERM_ALLOC
: Tests
{
63 my $eigen = gsl_eigen_genherm_alloc
(5);
64 isa_ok
($eigen, 'Math::GSL::Eigen');
67 sub GSL_EIGEN_GENHERMV_ALLOC
: Tests
{
68 my $eigen = gsl_eigen_genhermv_alloc
(5);
69 isa_ok
($eigen, 'Math::GSL::Eigen');
72 sub GSL_EIGEN_GEN_ALLOC
: Tests
{
73 my $eigen = gsl_eigen_gen_alloc
(5);
74 isa_ok
($eigen, 'Math::GSL::Eigen');
77 sub GSL_EIGEN_GENV_ALLOC
: Tests
{
78 my $eigen = gsl_eigen_genv_alloc
(5);
79 isa_ok
($eigen, 'Math::GSL::Eigen');
82 sub GSL_EIGEN_SYMM
: Tests
{
84 my $m = gsl_matrix_alloc
(2,2);
85 gsl_matrix_set_identity
($m);
86 my $v = gsl_vector_alloc
(2);
87 ok_status
(gsl_eigen_symm
($m, $v, $self->{eigen
}));
88 map { is
(gsl_vector_get
($v, $_), 1) } (0..1);
91 sub GSL_EIGEN_SYMMV
: Tests
{
92 my $w = gsl_eigen_symmv_alloc
(2);
93 my $m = gsl_matrix_alloc
(2,2);
94 gsl_matrix_set
($m, 0, 0, 2);
95 gsl_matrix_set
($m, 0, 1, 1);
96 gsl_matrix_set
($m, 1, 0, 1);
97 gsl_matrix_set
($m, 1, 1, 2);
98 my $eval = gsl_vector_alloc
(2);
99 my $evec = gsl_matrix_alloc
(2,2);
100 ok_status
(gsl_eigen_symmv
($m, $eval, $evec, $w));
101 is
(gsl_vector_get
($eval, 0), 3);
102 is
(gsl_vector_get
($eval, 1), 1);
103 my $x = gsl_matrix_get
($evec, 0, 0);
104 is
(gsl_matrix_get
($evec, 0, 1), -$x); #this is the eigenvector for the eigenvalue 1, which is the second eigenvalue in the $eval vector, but the GSL documentation says the first eigenvector should correspond to the first eigenvalue... where'e the error?
105 is
(sqrt($x**2+$x**2), 1);
107 $x = gsl_matrix_get
($evec, 1, 0);
108 is
(gsl_matrix_get
($evec, 1, 1), $x);
109 is
(sqrt($x**2+$x**2), 1);
111 my $v1 = gsl_vector_alloc
(2);
112 my $v2 = gsl_vector_alloc
(2);
113 gsl_matrix_get_col
($v1, $evec, 0);
114 gsl_matrix_get_col
($v2, $evec, 1);
115 gsl_vector_mul
($v1, $v2);
116 is
(gsl_vector_get
($v1, 0) + gsl_vector_get
($v1, 1) , 0);
119 sub GSL_EIGEN_SYMMV_SORT
: Tests
{
120 my $w = gsl_eigen_symmv_alloc
(2);
121 my $m = gsl_matrix_alloc
(2,2);
122 gsl_matrix_set
($m, 0, 0, 2);
123 gsl_matrix_set
($m, 0, 1, 1);
124 gsl_matrix_set
($m, 1, 0, 1);
125 gsl_matrix_set
($m, 1, 1, 2);
126 my $eval = gsl_vector_alloc
(2);
127 my $evec = gsl_matrix_alloc
(2,2);
128 ok_status
(gsl_eigen_symmv
($m, $eval, $evec, $w));
129 ok_status
(gsl_eigen_symmv_sort
($eval, $evec, $GSL_EIGEN_SORT_VAL_ASC));
130 is
(gsl_vector_get
($eval, 0), 1);
131 is
(gsl_vector_get
($eval, 1), 3);
132 my $x = gsl_matrix_get
($evec, 0, 0);
133 ok_similar
(gsl_matrix_get
($evec, 0, 1), -$x);
134 ok_similar
(sqrt($x**2+$x**2), 1);
136 $x = gsl_matrix_get
($evec, 1, 0);
137 ok_similar
(gsl_matrix_get
($evec, 1, 1), $x);
138 ok_similar
(sqrt($x**2+$x**2), 1);
140 my $v1 = gsl_vector_alloc
(2);
141 my $v2 = gsl_vector_alloc
(2);
142 gsl_matrix_get_col
($v1, $evec, 0);
143 gsl_matrix_get_col
($v2, $evec, 1);
144 gsl_vector_mul
($v1, $v2);
145 is
(gsl_vector_get
($v1, 0) + gsl_vector_get
($v1, 1) , 0);
148 sub GSL_EIGEN_HERM
: Tests
{
149 my $matrix = gsl_matrix_complex_alloc
(2, 2);
150 my $complex = gsl_complex_rect
(3,0);
151 gsl_matrix_complex_set
($matrix, 0, 0, $complex);
153 $complex = gsl_complex_rect
(2,1);
154 gsl_matrix_complex_set
($matrix, 0, 1, $complex);
156 $complex = gsl_complex_rect
(2,-1);
157 gsl_matrix_complex_set
($matrix, 1, 0, $complex);
159 $complex = gsl_complex_rect
(1,0);
160 gsl_matrix_complex_set
($matrix, 1, 1, $complex);
162 my $eigen = gsl_eigen_herm_alloc
(2);
163 my $vector = gsl_vector_alloc
(2);
164 is
(gsl_eigen_herm
($matrix, $vector, $eigen), 0);
165 is
(gsl_vector_get
($vector, 0), 2+sqrt(6));
166 is
(gsl_vector_get
($vector, 1), 2-sqrt(6));
168 sub GSL_EIGEN_HERMV
: Tests
{
169 my $matrix = gsl_matrix_complex_alloc
(2, 2);
170 my $complex = gsl_complex_rect
(3,0);
171 gsl_matrix_complex_set
($matrix, 0, 0, $complex);
173 $complex = gsl_complex_rect
(2,1);
174 gsl_matrix_complex_set
($matrix, 0, 1, $complex);
176 $complex = gsl_complex_rect
(2,-1);
177 gsl_matrix_complex_set
($matrix, 1, 0, $complex);
179 $complex = gsl_complex_rect
(1,0);
180 gsl_matrix_complex_set
($matrix, 1, 1, $complex);
182 my $eigen = gsl_eigen_hermv_alloc
(2);
183 my $vector = gsl_vector_alloc
(2);
184 my $evec = gsl_matrix_complex_alloc
(2,2);
185 is
(gsl_eigen_hermv
($matrix, $vector, $evec, $eigen), 0);
186 is
(gsl_vector_get
($vector, 0), 2+sqrt(6));
187 is
(gsl_vector_get
($vector, 1), 2-sqrt(6));
189 # my $x = gsl_matrix_complex_get($evec, 1, 1);
190 # ok_similar( [ gsl_parts($x) ], [-0.750532688285823, 0 ], "last row");
191 # ok_similar( [ gsl_parts($x) ], [ 2/(-1+sqrt(6)), 1/(-1+sqrt(6)) ], "first value");
195 sub GSL_EIGEN_NONSYMM
: Tests
{
196 my $matrix = gsl_matrix_alloc
(2, 2);
197 gsl_matrix_set
($matrix, 0, 0, -12);
198 gsl_matrix_set
($matrix, 1, 0, 7);
199 gsl_matrix_set
($matrix, 0, 1, 65);
200 gsl_matrix_set
($matrix, 1, 1, 59);
202 my $eigen = gsl_eigen_nonsymm_alloc
(2);
203 my $vector = gsl_vector_complex_alloc
(2);
204 is
(gsl_eigen_nonsymm
($matrix, $vector, $eigen), 0);
206 my $x = gsl_vector_complex_real
($vector);
207 my $y = gsl_vector_complex_imag
($vector);
209 # this interface seems hokey
210 is_similar
( gsl_vector_get
($x->{vector
}, 1), (47/2)+(0.5*sqrt(6861)) );
211 is_similar
( gsl_vector_get
($y->{vector
}, 1), 0 );
214 sub GSL_EIGEN_NONSYMM_Z
: Tests
{
215 my $matrix = gsl_matrix_alloc
(2, 2);
216 gsl_matrix_set
($matrix, 0, 0, -12);
217 gsl_matrix_set
($matrix, 1, 0, 7);
218 gsl_matrix_set
($matrix, 0, 1, 65);
219 gsl_matrix_set
($matrix, 1, 1, 59);
221 my $eigen = gsl_eigen_nonsymm_alloc
(2);
222 my $vector = gsl_vector_complex_alloc
(2);
223 my $Z = gsl_matrix_alloc
(2,2);
224 is
(gsl_eigen_nonsymm
($matrix, $vector, $eigen), 0);
225 is
(gsl_eigen_nonsymm_Z
($matrix,$vector, $Z, $eigen), 0);
226 ok
(is_similar
(gsl_matrix_get
($Z, 0, 0), 0.9958842418254068860784291, 0.005));
227 ok_similar
(gsl_matrix_get
($Z, 0, 1), 0.09063430301952179629793610, "Z matrix", 0.1);
228 ok_similar
(gsl_matrix_get
($Z, 1, 1), 0.9958842418254068860784291, "Z matrix", 0.005);
229 ok_similar
(gsl_matrix_get
($Z, 1, 0), 0.09063430301952179629793610, "Z matrix", 0.1);
231 my $x = gsl_vector_complex_real
($vector);
232 my $y = gsl_vector_complex_imag
($vector);
234 # this interface seems hokey
235 is_similar
( gsl_vector_get
($x->{vector
}, 1), (47/2)+(0.5*sqrt(6861)) );
236 is_similar
( gsl_vector_get
($y->{vector
}, 1), 0 );
239 sub GSL_EIGEN_NONSYMMV_Z
: Tests
{
240 my $matrix = gsl_matrix_alloc
(2, 2);
241 gsl_matrix_set
($matrix, 0, 0, -12);
242 gsl_matrix_set
($matrix, 1, 0, 7);
243 gsl_matrix_set
($matrix, 0, 1, 65);
244 gsl_matrix_set
($matrix, 1, 1, 59);
246 my $evec = gsl_matrix_complex_alloc
(2,2);
247 my $eigen = gsl_eigen_nonsymmv_alloc
(2);
248 my $vector = gsl_vector_complex_alloc
(2);
249 my $Z = gsl_matrix_alloc
(2,2);
250 ok_status
(gsl_eigen_nonsymmv_Z
($matrix,$vector, $evec, $Z, $eigen));
253 my $x = gsl_vector_complex_real
($vector);
254 my $y = gsl_vector_complex_imag
($vector);
256 # this interface seems hokey
257 is_similar
( gsl_vector_get
($x->{vector
}, 1), (47/2)+(0.5*sqrt(6861)) );
258 is_similar
( gsl_vector_get
($y->{vector
}, 1), 0 );
260 $x = gsl_matrix_complex_get
($evec, 1, 0);
261 ok_similar
(gsl_imag
($x), 0, "evec matrix");
262 ok_similar
(gsl_real
($x), 7/((71/2)+(.5*sqrt(6861))), "evec matrix", 0.01);
264 $x = gsl_matrix_complex_get
($evec, 0, 0);
265 ok_similar
(gsl_imag
($x), 0, "evec matrix");
266 ok_similar
(gsl_real
($x), 7/((71/2)-(.5*sqrt(6861))), "evec matrix", 0.19);
268 $x = gsl_matrix_complex_get
($evec, 0, 1);
269 $y = gsl_matrix_complex_get
($evec, 1, 1);
271 ok_similar
(gsl_imag
($x), 0, "evec matrix");
272 ok_similar
(gsl_imag
($y), 0, "evec matrix");
274 local $TODO = "matlab differences";
276 ok_similar
(gsl_real
($x), 1); # this is the value I get with maple
277 ok_similar
(gsl_real
($y), 1); # this is the value I get with maple
279 ok_similar
([ gsl_matrix_get
($Z, 0, 0)], [0.9958842418254068860784291] );
280 ok_similar
([ gsl_matrix_get
($Z, 0, 1)], [0.09063430301952179629793610] );
281 ok_similar
([ gsl_matrix_get
($Z, 1, 1)], [0.9958842418254068860784291] );
282 ok_similar
([ gsl_matrix_get
($Z, 1, 0)], [0.09063430301952179629793610] );