"Fossies" - the Fresh Open Source Software Archive  

Source code changes of the file "GENERATED/PDL/MatrixOps.pm" between
PDL-2.082.tar.gz and PDL-2.083.tar.gz

About: PDL (Perl Data Language) aims to turn perl into an efficient numerical language for scientific computing (similar to IDL and MatLab).

MatrixOps.pm  (PDL-2.082):MatrixOps.pm  (PDL-2.083)
skipping to change at line 17 skipping to change at line 17
our %EXPORT_TAGS = (Func=>\@EXPORT_OK); our %EXPORT_TAGS = (Func=>\@EXPORT_OK);
use PDL::Core; use PDL::Core;
use PDL::Exporter; use PDL::Exporter;
use DynaLoader; use DynaLoader;
our @ISA = ( 'PDL::Exporter','DynaLoader' ); our @ISA = ( 'PDL::Exporter','DynaLoader' );
push @PDL::Core::PP, __PACKAGE__; push @PDL::Core::PP, __PACKAGE__;
bootstrap PDL::MatrixOps ; bootstrap PDL::MatrixOps ;
#line 8 "matrixops.pd" #line 9 "matrixops.pd"
use strict; use strict;
use warnings; use warnings;
=head1 NAME =head1 NAME
PDL::MatrixOps -- Some Useful Matrix Operations PDL::MatrixOps -- Some Useful Matrix Operations
=head1 SYNOPSIS =head1 SYNOPSIS
skipping to change at line 131 skipping to change at line 131
=cut =cut
use Carp; use Carp;
use strict; use strict;
#line 134 "MatrixOps.pm" #line 134 "MatrixOps.pm"
=head1 FUNCTIONS =head1 FUNCTIONS
=cut =cut
#line 123 "matrixops.pd" #line 124 "matrixops.pd"
=head2 identity =head2 identity
=for sig =for sig
Signature: (n; [o]a(n,n)) Signature: (n; [o]a(n,n))
=for ref =for ref
Return an identity matrix of the specified size. If you hand in a Return an identity matrix of the specified size. If you hand in a
skipping to change at line 161 skipping to change at line 161
!(my $was_pdl = UNIVERSAL::isa($n,'PDL')) ? zeroes($n,$n) : !(my $was_pdl = UNIVERSAL::isa($n,'PDL')) ? zeroes($n,$n) :
$n->getndims == 0 ? zeroes($n->type, $n->at(0),$n->at(0)) : $n->getndims == 0 ? zeroes($n->type, $n->at(0),$n->at(0)) :
undef; undef;
if (!defined $out) { if (!defined $out) {
my @dims = $n->dims; my @dims = $n->dims;
$out = zeroes($n->type, @dims[0, 0, 2..$#dims]); $out = zeroes($n->type, @dims[0, 0, 2..$#dims]);
} }
(my $tmp = $out->diagonal(0,1))++; # work around perl -d "feature" (my $tmp = $out->diagonal(0,1))++; # work around perl -d "feature"
$was_pdl ? bless $out, ref($n) : $out; $was_pdl ? bless $out, ref($n) : $out;
} }
#line 178 "MatrixOps.pm"
#line 156 "matrixops.pd"
#line 157 "matrixops.pd"
=head2 stretcher =head2 stretcher
=for sig =for sig
Signature: (a(n); [o]b(n,n)) Signature: (a(n); [o]b(n,n))
=for usage =for usage
$mat = stretcher($eigenvalues); $mat = stretcher($eigenvalues);
skipping to change at line 188 skipping to change at line 186
=cut =cut
sub stretcher { sub stretcher {
my $in = shift; my $in = shift;
my $out = zeroes($in->dim(0),$in->dims); my $out = zeroes($in->dim(0),$in->dims);
my $tmp; # work around for perl -d "feature" my $tmp; # work around for perl -d "feature"
($tmp = $out->diagonal(0,1)) += $in; ($tmp = $out->diagonal(0,1)) += $in;
$out; $out;
} }
#line 208 "MatrixOps.pm"
#line 187 "matrixops.pd"
#line 188 "matrixops.pd"
=head2 inv =head2 inv
=for sig =for sig
Signature: (a(m,m); sv opt ) Signature: (a(m,m); sv opt )
=for usage =for usage
$a1 = inv($a, {$opt}); $a1 = inv($a, {$opt});
skipping to change at line 283 skipping to change at line 279
} }
my $out = lu_backsub($lu,$perm,$par,identity($x))->transpose->sever; my $out = lu_backsub($lu,$perm,$par,identity($x))->transpose->sever;
return $out return $out
unless($x->is_inplace); unless($x->is_inplace);
$x .= $out; $x .= $out;
$x; $x;
} }
#line 307 "MatrixOps.pm"
#line 287 "matrixops.pd"
#line 288 "matrixops.pd"
=head2 det =head2 det
=for sig =for sig
Signature: (a(m,m); sv opt) Signature: (a(m,m); sv opt)
=for usage =for usage
$det = det($a,{opt}); $det = det($a,{opt});
skipping to change at line 338 skipping to change at line 332
my($lu,$perm,$par); my($lu,$perm,$par);
if(exists ($opt->{lu}) and (ref $opt->{lu} eq 'ARRAY')) { if(exists ($opt->{lu}) and (ref $opt->{lu} eq 'ARRAY')) {
($lu,$perm,$par) = @{$opt->{lu}}; ($lu,$perm,$par) = @{$opt->{lu}};
} else { } else {
($lu,$perm,$par) = lu_decomp($x); ($lu,$perm,$par) = lu_decomp($x);
$opt->{lu} = [$lu,$perm,$par] $opt->{lu} = [$lu,$perm,$par]
if(exists($opt->{lu})); if(exists($opt->{lu}));
} }
defined $lu ? $lu->diagonal(0,1)->prodover * $par : PDL->zeroes(sbyte,1); defined $lu ? $lu->diagonal(0,1)->prodover * $par : PDL->zeroes(sbyte,1);
} }
#line 365 "MatrixOps.pm"
#line 346 "matrixops.pd"
#line 347 "matrixops.pd"
=head2 determinant =head2 determinant
=for sig =for sig
Signature: (a(m,m)) Signature: (a(m,m))
=for usage =for usage
$det = determinant($x); $det = determinant($x);
skipping to change at line 422 skipping to change at line 414
determinant($x->slice("0:".($i-1).",1:-1")-> determinant($x->slice("0:".($i-1).",1:-1")->
append($x->slice(($i+1).":-1,1:-1"))); append($x->slice(($i+1).":-1,1:-1")));
} }
# Do beginning and end submatrices # Do beginning and end submatrices
$sum += $x->slice("(0),(0)") * determinant($x->slice('1:-1,1:-1')); $sum += $x->slice("(0),(0)") * determinant($x->slice('1:-1,1:-1'));
$sum -= $x->slice("(-1),(0)") * determinant($x->slice('0:-2,1:-1')) * (1 - 2*( $n % 2)); $sum -= $x->slice("(-1),(0)") * determinant($x->slice('0:-2,1:-1')) * (1 - 2*( $n % 2));
return $sum; return $sum;
} }
#line 452 "MatrixOps.pm" #line 427 "MatrixOps.pm"
#line 958 "../../blib/lib/PDL/PP.pm"
=head2 eigens_sym =head2 eigens_sym
=for sig =for sig
Signature: ([phys]a(m); [o,phys]ev(n,n); [o,phys]e(n)) Signature: ([phys]a(m); [o,phys]ev(n,n); [o,phys]e(n))
=for ref =for ref
Eigenvalues and -vectors of a symmetric square matrix. If passed Eigenvalues and -vectors of a symmetric square matrix. If passed
skipping to change at line 468 skipping to change at line 458
($ev, $e) = eigens_sym($x); # e-vects & e-values ($ev, $e) = eigens_sym($x); # e-vects & e-values
$e = eigens_sym($x); # just eigenvalues $e = eigens_sym($x); # just eigenvalues
=for bad =for bad
eigens_sym ignores the bad-value flag of the input ndarrays. eigens_sym ignores the bad-value flag of the input ndarrays.
It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays.
=cut =cut
#line 506 "MatrixOps.pm"
#line 959 "../../blib/lib/PDL/PP.pm"
sub PDL::eigens_sym { sub PDL::eigens_sym {
my ($x) = @_; my ($x) = @_;
my (@d) = $x->dims; my (@d) = $x->dims;
barf "Need real square matrix for eigens_sym" barf "Need real square matrix for eigens_sym"
if $#d < 1 or $d[0] != $d[1]; if $#d < 1 or $d[0] != $d[1];
my ($n) = $d[0]; my ($n) = $d[0];
my ($sym) = 0.5*($x + $x->transpose); my ($sym) = 0.5*($x + $x->transpose);
my ($err) = PDL::max(abs($sym)); my ($err) = PDL::max(abs($sym));
barf "Need symmetric component non-zero for eigens_sym" barf "Need symmetric component non-zero for eigens_sym"
skipping to change at line 504 skipping to change at line 491
)->copy; )->copy;
my $ev = PDL->zeroes($sym->dims); my $ev = PDL->zeroes($sym->dims);
my $e = PDL->zeroes($sym->index(0)->dims); my $e = PDL->zeroes($sym->index(0)->dims);
&PDL::_eigens_sym_int($lt, $ev, $e); &PDL::_eigens_sym_int($lt, $ev, $e);
return $ev->transpose, $e return $ev->transpose, $e
if(wantarray); if(wantarray);
$e; #just eigenvalues $e; #just eigenvalues
} }
#line 545 "MatrixOps.pm"
#line 960 "../../blib/lib/PDL/PP.pm"
*eigens_sym = \&PDL::eigens_sym; *eigens_sym = \&PDL::eigens_sym;
#line 552 "MatrixOps.pm"
#line 958 "../../blib/lib/PDL/PP.pm"
=head2 eigens =head2 eigens
=for sig =for sig
Signature: ([phys]a(m); [o,phys]ev(l,n,n); [o,phys]e(l,n)) Signature: ([phys]a(m); [o,phys]ev(l,n,n); [o,phys]e(l,n))
=for ref =for ref
Real eigenvalues and -vectors of a real square matrix. Real eigenvalues and -vectors of a real square matrix.
skipping to change at line 581 skipping to change at line 562
($ev, $e) = eigens($x); # e'vects & e'vals ($ev, $e) = eigens($x); # e'vects & e'vals
$e = eigens($x); # just eigenvalues $e = eigens($x); # just eigenvalues
=for bad =for bad
eigens ignores the bad-value flag of the input ndarrays. eigens ignores the bad-value flag of the input ndarrays.
It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays.
=cut =cut
#line 631 "MatrixOps.pm"
#line 959 "../../blib/lib/PDL/PP.pm"
sub PDL::eigens { sub PDL::eigens {
my ($x) = @_; my ($x) = @_;
my (@d) = $x->dims; my (@d) = $x->dims;
my $n = $d[0]; my $n = $d[0];
barf "Need real square matrix for eigens" barf "Need real square matrix for eigens"
if $#d < 1 or $d[0] != $d[1]; if $#d < 1 or $d[0] != $d[1];
my $deviation = PDL::max(abs($x - $x->transpose))/PDL::max(abs($x)); my $deviation = PDL::max(abs($x - $x->transpose))/PDL::max(abs($x));
if ( $deviation <= 1e-5 ) { if ( $deviation <= 1e-5 ) {
#taken from eigens_sym code #taken from eigens_sym code
skipping to change at line 630 skipping to change at line 608
my $ev = PDL->zeroes(2, $x->dims); my $ev = PDL->zeroes(2, $x->dims);
my $e = PDL->zeroes(2, $x->index(0)->dims); my $e = PDL->zeroes(2, $x->index(0)->dims);
&PDL::_eigens_int($x->clump(0,1), $ev, $e); &PDL::_eigens_int($x->clump(0,1), $ev, $e);
return $ev->index(0)->transpose->sever, $e->index(0)->sever return $ev->index(0)->transpose->sever, $e->index(0)->sever
if(wantarray); if(wantarray);
return $e->index(0)->sever; #just eigenvalues return $e->index(0)->sever; #just eigenvalues
} }
} }
#line 683 "MatrixOps.pm"
#line 960 "../../blib/lib/PDL/PP.pm"
*eigens = \&PDL::eigens; *eigens = \&PDL::eigens;
#line 690 "MatrixOps.pm"
#line 958 "../../blib/lib/PDL/PP.pm"
=head2 svd =head2 svd
=for sig =for sig
Signature: (a(n,m); [t]w(wsize); [o]u(n,m); [o,phys]z(n); [o]v(n,n)) Signature: (a(n,m); [t]w(wsize); [o]u(n,m); [o,phys]z(n); [o]v(n,n))
=for usage =for usage
($u, $s, $v) = svd($x); ($u, $s, $v) = svd($x);
skipping to change at line 699 skipping to change at line 671
$r2 *= $s; # implicit broadcasting for cheap mult. $r2 *= $s; # implicit broadcasting for cheap mult.
$x .= $r2 x $r1; # a gets r2 x ess x r1 $x .= $r2 x $r1; # a gets r2 x ess x r1
} }
=for bad =for bad
svd ignores the bad-value flag of the input ndarrays. svd ignores the bad-value flag of the input ndarrays.
It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays.
=cut =cut
#line 761 "MatrixOps.pm"
#line 960 "../../blib/lib/PDL/PP.pm"
*svd = \&PDL::svd; *svd = \&PDL::svd;
#line 768 "MatrixOps.pm"
#line 800 "matrixops.pd" #line 801 "matrixops.pd"
=head2 lu_decomp =head2 lu_decomp
=for sig =for sig
Signature: (a(m,m); [o]lu(m,m); [o]perm(m); [o]parity) Signature: (a(m,m); [o]lu(m,m); [o]perm(m); [o]parity)
=for ref =for ref
LU decompose a matrix, with row permutation LU decompose a matrix, with row permutation
skipping to change at line 880 skipping to change at line 848
return if !$notbig->isempty; return if !$notbig->isempty;
# Divide by the diagonal element (which is now the largest element) # Divide by the diagonal element (which is now the largest element)
my $tout; my $tout;
($tout = $out->slice("($col),".($col+1).":$n1")) /= $big->slice('*1'); ($tout = $out->slice("($col),".($col+1).":$n1")) /= $big->slice('*1');
} # end of pivoting part } # end of pivoting part
} # end of column loop } # end of column loop
wantarray ? ($out,$permute,$parity) : $out; wantarray ? ($out,$permute,$parity) : $out;
} }
#line 947 "MatrixOps.pm"
#line 979 "matrixops.pd"
#line 980 "matrixops.pd"
=head2 lu_decomp2 =head2 lu_decomp2
=for sig =for sig
Signature: (a(m,m); [o]lu(m,m)) Signature: (a(m,m); [o]lu(m,m))
=for ref =for ref
LU decompose a matrix, with no row permutation LU decompose a matrix, with no row permutation
skipping to change at line 998 skipping to change at line 964
if($col < $n1) { if($col < $n1) {
# Divide the rest of the column by the diagonal element # Divide the rest of the column by the diagonal element
my $tmp; # work around for perl -d "feature" my $tmp; # work around for perl -d "feature"
($tmp = $out->slice("($col),".($col+1).":$n1")) /= $diagonal->index($col)- >dummy(0,$n1-$col); ($tmp = $out->slice("($col),".($col+1).":$n1")) /= $diagonal->index($col)- >dummy(0,$n1-$col);
} }
} # end of column loop } # end of column loop
wantarray ? ($out,$perm,$par) : $out; wantarray ? ($out,$perm,$par) : $out;
} }
#line 1069 "MatrixOps.pm"
#line 1101 "matrixops.pd"
#line 1102 "matrixops.pd"
=head2 lu_backsub =head2 lu_backsub
=for sig =for sig
Signature: (lu(m,m); perm(m); b(m)) Signature: (lu(m,m); perm(m); b(m))
=for ref =for ref
Solve A x = B for matrix A, by back substitution into A's LU decomposition. Solve A x = B for matrix A, by back substitution into A's LU decomposition.
skipping to change at line 1213 skipping to change at line 1177
)->sumover; )->sumover;
($tmp = $out->index($r1)) /= $ludiag->index($r1)->dummy(0); # TODO: check broadcast dims ($tmp = $out->index($r1)) /= $ludiag->index($r1)->dummy(0); # TODO: check broadcast dims
} }
if ($y->is_inplace) { if ($y->is_inplace) {
$y->setdims([$out->dims]) if !PDL::all($y->shape == $out->shape); # assgn n eeds same shape $y->setdims([$out->dims]) if !PDL::all($y->shape == $out->shape); # assgn n eeds same shape
$y .= $out; $y .= $out;
} }
$out; $out;
} }
#line 1287 "MatrixOps.pm" #line 1219 "MatrixOps.pm"
#line 958 "../../blib/lib/PDL/PP.pm"
=head2 simq =head2 simq
=for sig =for sig
Signature: ([phys]a(n,n); [phys]b(n); [o,phys]x(n); int [o,phys]ips(n); int fl ag) Signature: ([phys]a(n,n); [phys]b(n); [o,phys]x(n); int [o,phys]ips(n); int fl ag)
=for ref =for ref
Solution of simultaneous linear equations, C<a x = b>. Solution of simultaneous linear equations, C<a x = b>.
skipping to change at line 1251 skipping to change at line 1213
See also L</lu_backsub>, which does the same thing with a slightly See also L</lu_backsub>, which does the same thing with a slightly
less opaque interface. less opaque interface.
=for bad =for bad
simq ignores the bad-value flag of the input ndarrays. simq ignores the bad-value flag of the input ndarrays.
It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays.
=cut =cut
#line 1333 "MatrixOps.pm"
#line 960 "../../blib/lib/PDL/PP.pm"
*simq = \&PDL::simq; *simq = \&PDL::simq;
#line 1340 "MatrixOps.pm"
#line 958 "../../blib/lib/PDL/PP.pm"
=head2 squaretotri =head2 squaretotri
=for sig =for sig
Signature: (a(n,n); b(m)) Signature: (a(n,n); [o]b(m))
=for ref =for ref
Convert a symmetric square matrix to triangular vector storage. Convert a lower-triangular square matrix to triangular vector storage.
Ignores upper half of input.
=for bad =for bad
squaretotri does not process bad values. squaretotri does not process bad values.
It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays.
=cut =cut
#line 1368 "MatrixOps.pm"
#line 960 "../../blib/lib/PDL/PP.pm"
*squaretotri = \&PDL::squaretotri; *squaretotri = \&PDL::squaretotri;
#line 1375 "MatrixOps.pm"
#line 1396 "matrixops.pd" #line 1390 "matrixops.pd"
=head1 AUTHOR =head1 AUTHOR
Copyright (C) 2002 Craig DeForest (deforest@boulder.swri.edu), Copyright (C) 2002 Craig DeForest (deforest@boulder.swri.edu),
R.J.R. Williams (rjrw@ast.leeds.ac.uk), Karl Glazebrook R.J.R. Williams (rjrw@ast.leeds.ac.uk), Karl Glazebrook
(kgb@aaoepp.aao.gov.au). There is no warranty. You are allowed to (kgb@aaoepp.aao.gov.au). There is no warranty. You are allowed to
redistribute and/or modify this work under the same conditions as PDL redistribute and/or modify this work under the same conditions as PDL
itself. If this file is separated from the PDL distribution, then the itself. If this file is separated from the PDL distribution, then the
PDL copyright notice should be included in this file. PDL copyright notice should be included in this file.
=cut =cut
#line 1394 "MatrixOps.pm" #line 1307 "MatrixOps.pm"
# Exit with OK status # Exit with OK status
1; 1;
 End of changes. 33 change blocks. 
63 lines changed or deleted 16 lines changed or added

Home  |  About  |  Features  |  All  |  Newest  |  Dox  |  Diffs  |  RSS Feeds  |  Screenshots  |  Comments  |  Imprint  |  Privacy  |  HTTP(S)