Source

perl-PDL-Lib-Linear-Solve / matrix_inversion.pl

#!/usr/bin/perl -w

use strict;

use MatInv;
use PDL;

my ($a, $v);
my ($new_a, $new_v);
sub mytest
{
    print "\n\n-------------\n\n";
print "\$a = " , $a;
print "\$v = ", $v, "\n";    
($new_a, $new_v) = linear_solver($a,$v);
print "\$new_a = ", $new_a;
print "\$new_v = ", $new_v;
print "\n\n";
    
}

sub eye
{
    my $n = shift;
    my $pdl = zeroes($n,$n);
    $pdl->diagonal(0,1)++;
    return $pdl;
}

$a = pdl([[0,2,3], [4, 5 ,6], [7, 8, 9]]);
$v = pdl([1, 5 ,10]);

mytest();

$a = pdl([[1,0,1],[0,1,0],[0,0,1]]);
$v = pdl([3,2,1]);

mytest();

$a = pdl([[1,0,1],[0,1,0],[0,0,1]]);
$v = eye(3);
mytest();
print $new_a*$new_v;

$a = pdl([[1,0,1],[0,1,0],[0,0,1]]);
$v = pdl([
         [[1,2,3],[4,5,6],[7,8,9]],
         [[10,20,30],[40,50,60],[70,80,90]],
         [[100,200,300],[400,500,600],[700,800,900]]
         ]);
mytest();