[system] / trunk / pg / lib / Value / Matrix.pm Repository:
ViewVC logotype

View of /trunk/pg/lib/Value/Matrix.pm

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3171 - (download) (as text) (annotate)
Tue Feb 15 21:53:23 2005 UTC (15 years ago) by dpvc
File size: 11659 byte(s)
Improved the Real(), Complex(), Point(), Vector(), Matrix() and
String() constructors so that they will process formulas passed to
them as strings rather than requiring perl objects for these.

For example, you can use Real("2/3") rather than Real(2/3) if you
want.  Also, Real("1+x") will return a formula returning a real
(essentially the same as Formula("1+x") in this case).

    1 ###########################################################################
    2 #
    3 #  Implements the Matrix class.
    4 #
    5 #    @@@ Still needs lots of work @@@
    6 #
    7 package Value::Matrix;
    8 my $pkg = 'Value::Matrix';
    9 
   10 use strict;
   11 use vars qw(@ISA);
   12 @ISA = qw(Value);
   13 
   14 use overload
   15        '+'   => \&add,
   16        '-'   => \&sub,
   17        '*'   => \&mult,
   18        '/'   => \&div,
   19        '**'  => \&power,
   20        '.'   => \&Value::_dot,
   21        'x'   => \&Value::cross,
   22        '<=>' => \&compare,
   23        'cmp' => \&Value::cmp,
   24        'neg' => sub {$_[0]->neg},
   25   'nomethod' => \&Value::nomethod,
   26         '""' => \&stringify;
   27 
   28 #
   29 #  Convert a value to a matrix.  The value can be:
   30 #     a list of numbers or list of (nested) references to arrays of numbers
   31 #     a point, vector or matrix object, a matrix-valued formula, or a string
   32 #     that evaluates to a matrix
   33 #
   34 sub new {
   35   my $self = shift; my $class = ref($self) || $self;
   36   my $M = shift; $M = Value::makeValue($M) if !ref($M) && scalar(@_) == 0;
   37   return bless {data => $M->data}, $class
   38     if (Value::class($M) =~ m/Point|Vector|Matrix/ && scalar(@_) == 0);
   39   return $M if (Value::isFormula($M) && $M->type eq Value::class($self));
   40   $M = [$M,@_] if (ref($M) ne 'ARRAY' || scalar(@_) > 0);
   41   Value::Error("Matrices must have at least one entry") unless scalar(@{$M}) > 0;
   42   return $self->matrixMatrix(@{$M}) if ref($M->[0]) =~ m/ARRAY|Matrix/;
   43   return $self->numberMatrix(@{$M});
   44 }
   45 
   46 #
   47 #  (Recursively) make a matrix from a list of array refs
   48 #  and report errors about the entry types
   49 #
   50 sub matrixMatrix {
   51   my $self = shift; my $class = ref($self) || $self;
   52   my ($x,$m); my @M = (); my $isFormula = 0;
   53   foreach $x (@_) {
   54     if (Value::isFormula($x)) {push(@M,$x); $isFormula = 1} else {
   55       $m = $pkg->new($x); push(@M,$m);
   56       $isFormula = 1 if Value::isFormula($m);
   57     }
   58   }
   59   my ($type,$len) = ($M[0]->entryType->{name},$M[0]->length);
   60   foreach $x (@M) {
   61     Value::Error("Matrix rows must all be the same type")
   62       unless (defined($x->entryType) && $type eq $x->entryType->{name});
   63     Value::Error("Matrix rows must all be the same length") unless ($len eq $x->length);
   64   }
   65   return $self->formula([@M]) if $isFormula;
   66   bless {data => [@M]}, $class;
   67 }
   68 
   69 #
   70 #  Form a 1 x n matrix from a list of numbers
   71 #  (could become a row of an  m x n  matrix)
   72 #
   73 sub numberMatrix {
   74   my $self = shift; my $class = ref($self) || $self;
   75   my @M = (); my $isFormula = 0;
   76   foreach my $x (@_) {
   77     $x = Value::makeValue($x);
   78     Value::Error("Matrix row entries must be numbers") unless Value::isNumber($x);
   79     push(@M,$x); $isFormula = 1 if Value::isFormula($x);
   80   }
   81   return $self->formula([@M]) if $isFormula;
   82   bless {data => [@M]}, $class;
   83 }
   84 
   85 #
   86 #  Recursively get the entries in the matrix and return
   87 #  an array of (references to arrays of ... ) numbers
   88 #
   89 sub value {
   90   my $self = shift;
   91   my $M = $self->data;
   92   return @{$M} if Value::class($M->[0]) ne 'Matrix';
   93   my @M = ();
   94   foreach my $x (@{$M}) {push(@M,[$x->value])}
   95   return @M;
   96 }
   97 #
   98 #  The number of rows in the matrix (for n x m)
   99 #  or the number of entries in a 1 x n matrix
  100 #
  101 sub length {return scalar(@{shift->{data}})}
  102 #
  103 #  Recursively get the dimensions of the matrix.
  104 #  Returns (n) for a 1 x n, or (n,m) for an n x m, etc.
  105 #
  106 sub dimensions {
  107   my $self = shift;
  108   my $r = $self->length;
  109   my $v = $self->data;
  110   return ($r,) if (Value::class($v->[0]) ne 'Matrix');
  111   return ($r,$v->[0]->dimensions);
  112 }
  113 #
  114 #  Return the proper type for the matrix
  115 #
  116 sub typeRef {
  117   my $self = shift;
  118   return Value::Type($self->class, $self->length, $Value::Type{number})
  119     if (Value::class($self->data->[0]) ne 'Matrix');
  120   return Value::Type($self->class, $self->length, $self->data->[0]->typeRef);
  121 }
  122 
  123 #
  124 #  True if the matrix is a square matrix
  125 #
  126 sub isSquare {
  127   my $self = shift;
  128   my @d = $self->dimensions;
  129   return 0 if scalar(@d) > 2;
  130   return 1 if scalar(@d) == 1 && $d[0] == 1;
  131   return $d[0] == $d[1];
  132 }
  133 
  134 #
  135 #  True if the matrix is 1-dimensional (i.e., is a matrix row)
  136 #
  137 sub isRow {
  138   my $self = shift;
  139   my @d = $self->dimensions;
  140   return scalar(@d) == 1;
  141 }
  142 
  143 #
  144 #  See if the matrix is an Indenity matrix
  145 #
  146 sub isOne {
  147   my $self = shift;
  148   return 0 unless $self->isSquare;
  149   my $i = 0;
  150   foreach my $row (@{$self->{data}}) {
  151     my $j = 0;
  152     foreach my $k (@{$row->{data}}) {
  153       return 0 unless $k eq (($i == $j)? "1": "0");
  154       $j++;
  155     }
  156     $i++;
  157   }
  158   return 1;
  159 }
  160 
  161 #
  162 #  See if the matrix is all zeros
  163 #
  164 sub isZero {
  165   my $self = shift;
  166   foreach my $x (@{$self->{data}}) {return 0 unless $x->isZero}
  167   return 1;
  168 }
  169 
  170 #
  171 #  Make arbitrary data into a matrix, if possible
  172 #
  173 sub promote {
  174   my $x = shift;
  175   return $pkg->new($x,@_) if scalar(@_) > 0 || ref($x) eq 'ARRAY';
  176   return $x if ref($x) eq $pkg;
  177   return $pkg->make(@{$x->data}) if Value::class($x) =~ m/Point|Vector/;
  178   Value::Error("Can't convert ".Value::showClass($x)." to a Matrix");
  179 }
  180 
  181 ############################################
  182 #
  183 #  Operations on matrices
  184 #
  185 
  186 sub add {
  187   my ($l,$r,$flag) = @_;
  188   if ($l->promotePrecedence($r)) {return $r->add($l,!$flag)}
  189   ($l,$r) = (promote($l)->data,promote($r)->data);
  190   Value::Error("Matrix addition with different dimensions")
  191     unless scalar(@{$l}) == scalar(@{$r});
  192   my @s = ();
  193   foreach my $i (0..scalar(@{$l})-1) {push(@s,$l->[$i] + $r->[$i])}
  194   return $pkg->make(@s);
  195 }
  196 
  197 sub sub {
  198   my ($l,$r,$flag) = @_;
  199   if ($l->promotePrecedence($r)) {return $r->sub($l,!$flag)}
  200   ($l,$r) = (promote($l)->data,promote($r)->data);
  201   Value::Error("Matrix subtraction with different dimensions")
  202     unless scalar(@{$l}) == scalar(@{$r});
  203   if ($flag) {my $tmp = $l; $l = $r; $r = $tmp};
  204   my @s = ();
  205   foreach my $i (0..scalar(@{$l})-1) {push(@s,$l->[$i] - $r->[$i])}
  206   return $pkg->make(@s);
  207 }
  208 
  209 sub mult {
  210   my ($l,$r,$flag) = @_;
  211   if ($l->promotePrecedence($r)) {return $r->mult($l,!$flag)}
  212   #
  213   #  Constant multiplication
  214   #
  215   if (Value::matchNumber($r) || Value::isComplex($r)) {
  216     my @coords = ();
  217     foreach my $x (@{$l->data}) {push(@coords,$x*$r)}
  218     return $pkg->make(@coords);
  219   }
  220   #
  221   #  Make points and vectors into columns if they are on the right
  222   #
  223   if (!$flag && Value::class($r) =~ m/Point|Vector/)
  224     {$r = (promote($r))->transpose} else {$r = promote($r)}
  225   #
  226   if ($flag) {my $tmp = $l; $l = $r; $r = $tmp}
  227   my @dl = $l->dimensions; my @dr = $r->dimensions;
  228   if (scalar(@dl) == 1) {@dl = (1,@dl); $l = $pkg->make($l)}
  229   if (scalar(@dr) == 1) {@dr = (@dr,1); $r = $pkg->make($r)->transpose}
  230   Value::Error("Can only multiply 2-dimensional matrices") if scalar(@dl) > 2 || scalar(@dr) > 2;
  231   Value::Error("Matices of dimensions $dl[0]x$dl[1] and $dr[0]x$dr[1] can't be multiplied")
  232     unless ($dl[1] == $dr[0]);
  233   #
  234   #  Do matrix multiplication
  235   #
  236   my @l = $l->value; my @r = $r->value;
  237   my @M = ();
  238   foreach my $i (0..$dl[0]-1) {
  239     my @row = ();
  240     foreach my $j (0..$dr[1]-1) {
  241       my $s = 0;
  242       foreach my $k (0..$dl[1]-1) {$s += $l[$i]->[$k] * $r[$k]->[$j]}
  243       push(@row,$s);
  244     }
  245     push(@M,$pkg->make(@row));
  246   }
  247   return $pkg->make(@M);
  248 }
  249 
  250 sub div {
  251   my ($l,$r,$flag) = @_;
  252   if ($l->promotePrecedence($r)) {return $r->div($l,!$flag)}
  253   Value::Error("Can't divide by a Matrix") if $flag;
  254   Value::Error("Matrices can only be divided by numbers")
  255     unless (Value::matchNumber($r) || Value::isComplex($r));
  256   Value::Error("Division by zero") if $r == 0;
  257   my @coords = ();
  258   foreach my $x (@{$l->data}) {push(@coords,$x/$r)}
  259   return $pkg->make(@coords);
  260 }
  261 
  262 sub power {
  263   my ($l,$r,$flag) = @_;
  264   if ($l->promotePrecedence($r)) {return $r->power($l,!$flag)}
  265   Value::Error("Can't use Matrices in exponents") if $flag;
  266   Value::Error("Only square matrices can be raised to a power") unless $l->isSquare;
  267   return Value::Matrix::I($l->length) if $r == 0;
  268   Value::Error("Matrix powers must be positive integers") unless $r =~ m/^[1-9]\d*$/;
  269   my $M = $l; foreach my $i (2..$r) {$M = $M*$l}
  270   return $M;
  271 }
  272 
  273 #
  274 #  Do lexicographic comparison
  275 #
  276 sub compare {
  277   my ($l,$r,$flag) = @_;
  278   if ($l->promotePrecedence($r)) {return $r->compare($l,!$flag)}
  279   ($l,$r) = (promote($l)->data,promote($r)->data);
  280   Value::Error("Matrix comparison with different dimensions")
  281     unless scalar(@{$l}) == scalar(@{$r});
  282   if ($flag) {my $tmp = $l; $l = $r; $r = $tmp};
  283   my $cmp = 0;
  284   foreach my $i (0..scalar(@{$l})-1) {
  285     $cmp = $l->[$i] <=> $r->[$i];
  286     last if $cmp;
  287   }
  288   return $cmp;
  289 }
  290 
  291 sub neg {
  292   my $p = promote(@_)->data;
  293   my @coords = ();
  294   foreach my $x (@{$p}) {push(@coords,-$x)}
  295   return $pkg->make(@coords);
  296 }
  297 
  298 #
  299 #  Transpose an  n x m  matrix
  300 #
  301 sub transpose {
  302   my $self = shift;
  303   my @d = $self->dimensions;
  304   if (scalar(@d) == 1) {@d = (1,@d); $self = $pkg->make($self)}
  305   Value::Error("Can't transpose ".scalar(@d)."-dimensional matrices") unless scalar(@d) == 2;
  306   my @M = (); my $M = $self->data;
  307   foreach my $j (0..$d[1]-1) {
  308     my @row = ();
  309     foreach my $i (0..$d[0]-1) {push(@row,$M->[$i]->data->[$j])}
  310     push(@M,$pkg->make(@row));
  311   }
  312   return $pkg->make(@M);
  313 }
  314 
  315 #
  316 #  Get an identity matrix of the requested size
  317 #
  318 sub I {
  319   my $d = shift; $d = shift if ref($d) eq $pkg;
  320   my @M = (); my @Z = split('',0 x $d);
  321   foreach my $i (0..$d-1) {
  322     my @row = @Z; $row[$i] = 1;
  323     push(@M,$pkg->make(@row));
  324   }
  325   return $pkg->make(@M);
  326 }
  327 
  328 #
  329 #  Extract a given row from the matrix
  330 #
  331 sub row {
  332   my $M = promote(shift); my $i = shift;
  333   return if $i == 0; $i-- if $i > 0;
  334   if ($M->isRow) {return if $i != 0; return $M}
  335   return $M->data->[$i];
  336 }
  337 
  338 #
  339 #  Extract a given element from the matrix
  340 #
  341 sub element {
  342   my $M = promote(shift);
  343   return $M->extract(@_);
  344 }
  345 
  346 #
  347 #  Extract a given column from the matrix
  348 #
  349 sub column {
  350   my $M = promote(shift); my $j = shift;
  351   return if $j == 0; $j-- if $j > 0;
  352   my @d = $M->dimensions; my @col = ();
  353   return if $j+1 > $d[1];
  354   return $M->data->[$j] if scalar(@d) == 1;
  355   foreach my $row (@{$M->data}) {push(@col,$pkg->make($row->data->[$j]))}
  356   return $pkg->make(@col);
  357 }
  358 
  359 # @@@ removeRow, removeColumn @@@
  360 # @@@ Det, inverse @@@
  361 
  362 ############################################
  363 #
  364 #  Generate the various output formats
  365 #
  366 
  367 sub stringify {
  368   my $self = shift;
  369   return $self->TeX(undef,$self->{open},$self->{close}) if $$Value::context->flag('StringifyAsTeX');
  370   return $self->string(undef,$self->{open},$self->{close})
  371 }
  372 
  373 sub string {
  374   my $self = shift; my $equation = shift;
  375   my $def = ($equation->{context} || $$Value::context)->lists->get('Matrix');
  376   my $open = shift || $def->{open}; my $close = shift || $def->{close};
  377   my @coords = ();
  378   foreach my $x (@{$self->data}) {
  379     if (Value::isValue($x)) {push(@coords,$x->string($equation,$open,$close))}
  380       else {push(@coords,$x)}
  381   }
  382   return $open.join(',',@coords).$close;
  383 }
  384 
  385 #
  386 #  Use \matrix to lay out matrices
  387 #
  388 sub TeX {
  389   my $self = shift; my $equation = shift;
  390   my $def = ($equation->{context} || $$Value::context)->lists->get('Matrix');
  391   my $open = shift || $def->{open}; my $close = shift || $def->{close};
  392   $open = '\{' if $open eq '{'; $close = '\}' if $close eq '}';
  393   my $TeX = ''; my @entries = (); my $d;
  394   if ($self->isRow) {
  395     foreach my $x (@{$self->data}) {
  396       push(@entries,(Value::isValue($x))? $x->TeX($equation): $x);
  397     }
  398     $TeX .= join(' &',@entries) . "\n";
  399     $d = scalar(@entries);
  400   } else {
  401     foreach my $row (@{$self->data}) {
  402       foreach my $x (@{$row->data}) {
  403         push(@entries,(Value::isValue($x))? $x->TeX($equation): $x);
  404       }
  405       $TeX .= join(' &',@entries) . '\cr'."\n";
  406       $d = scalar(@entries); @entries = ();
  407     }
  408   }
  409   return '\left'.$open.'\begin{array}{'.('c'x$d).'}'."\n".$TeX.'\end{array}\right'.$close;
  410 }
  411 
  412 ###########################################################################
  413 
  414 1;
  415 

aubreyja at gmail dot com
ViewVC Help
Powered by ViewVC 1.0.9