1#!/usr/bin/perl 2######################### We start with some black magic to print on failure. 3use lib '../blib/lib','../blib/arch'; 4use strict; 5use warnings; 6use vars qw($loaded); 7 8BEGIN {$| = 1; print "1..74\n";} 9END {print "not ok 1\n" unless $loaded;} 10use Math::Cephes::Matrix qw(mat); 11$loaded = 1; 12print "ok 1\n"; 13 14######################### End of black magic. 15 16# util 17my $count = 1; 18my $eps = 1e-07; 19sub ok { 20 local($^W) = 0; 21 $count++; 22 my ($package, $file, $line) = caller; 23 my ($value, $true, $skip) = @_; 24 $skip ||= ''; 25 $skip = "# skip ($skip)" if $skip; 26 my $error = sprintf( "%12.8f", abs($value - $true)); 27 print($error < $eps ? "ok $count $skip\n" : 28 "not ok $count (expected $true: got $value) at $file line $line\n"); 29} 30 31my $M = Math::Cephes::Matrix->new([ [1, 2, -1], [2, -3, 1], [1, 0, 3]]); 32my $B = [2, -1, 10]; 33my $X = $M->simq($B); 34ok( $X->[0], 1); 35ok( $X->[1], 2); 36ok( $X->[2], 3); 37my $C = Math::Cephes::Matrix->new([ [1, 2, 4], [2, 9, 2], [6, 2, 7]]); 38my $I = $C->inv(); 39my $T = $I->mul($C)->coef; 40ok( $T->[0]->[0], 1); 41ok( $T->[1]->[1], 1); 42ok( $T->[2]->[2], 1); 43ok( $T->[0]->[1], 0); 44ok( $T->[1]->[0], 0); 45ok( $T->[2]->[0], 0); 46my $V = $M->mul($X); 47ok( $V->[0], $B->[0]); 48ok( $V->[1], $B->[1]); 49ok( $V->[2], $B->[2]); 50my $D = $M->add($C)->coef; 51ok( $D->[0]->[0], 2); 52ok( $D->[1]->[1], 6); 53ok( $D->[2]->[2], 10); 54ok( $D->[0]->[1], 4); 55ok( $D->[1]->[0], 4); 56ok( $D->[2]->[0], 7); 57$D = $M->sub($C)->coef; 58ok( $D->[0]->[0], 0); 59ok( $D->[1]->[1], -12); 60ok( $D->[2]->[2], -4); 61ok( $D->[0]->[1], 0); 62ok( $D->[1]->[0], 0); 63ok( $D->[2]->[0], -5); 64my $H = $C->transp()->coef; 65ok( $H->[0]->[0], 1); 66ok( $H->[1]->[1], 9); 67ok( $H->[2]->[2], 7); 68ok( $H->[0]->[1], 2); 69ok( $H->[1]->[0], 2); 70ok( $H->[2]->[0], 4); 71my $R = $M->div($C); 72my $Q = $R->mul($C)->coef; 73my $Mc = $M->coef; 74for (my $i=0; $i<3; $i++) { 75 for (my $j=0; $j<3; $j++) { 76 ok($Q->[$i]->[$j], $Mc->[$i]->[$j]); 77 } 78} 79$R = $M->mul($C)->coef; 80ok( $R->[0]->[0], -1); 81ok( $R->[1]->[1], -21); 82ok( $R->[2]->[2], 25); 83ok( $R->[0]->[1], 18); 84ok( $R->[1]->[0], 2); 85ok( $R->[2]->[0], 19); 86 87$C->clr(); 88$R = $C->coef; 89ok( $R->[0]->[0], 0); 90ok( $R->[2]->[2], 0); 91ok( $R->[1]->[0], 0); 92ok( $R->[2]->[0], 0); 93 94$C->clr(3); 95$R = $C->coef; 96ok( $R->[0]->[0], 3); 97ok( $R->[2]->[2], 3); 98ok( $R->[1]->[0], 3); 99ok( $R->[2]->[0], 3); 100 101my $S = Math::Cephes::Matrix->new([ [1, 2, 3], [2, 2, 3], [3, 3, 4]]); 102my ($E, $EV1) = $S->eigens(); 103my $EV = $EV1->coef; 104for (my $i=0; $i<3; $i++) { 105 my $v = []; 106 for (my $j=0; $j<3; $j++) { 107 $v->[$j] = $EV->[$i]->[$j]; 108 } 109 my $sv = $S->mul($v); 110 for (my $j=0; $j<3; $j++) { 111 ok($sv->[$j], $E->[$i]*$v->[$j]); 112 } 113} 114 115my $Z = $M->new()->coef; 116for (my $i=0; $i<3; $i++) { 117 for (my $j=0; $j<3; $j++) { 118 ok($Z->[$i]->[$j], $Mc->[$i]->[$j]); 119 } 120} 121$Z->[0]->[0] = 5; 122ok($Mc->[0]->[0], 1); 123ok($Z->[0]->[0], 5); 124 125