[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 2625 - (download) (as text) (annotate)
Mon Aug 16 18:35:12 2004 UTC (15 years, 3 months ago) by dpvc
File size: 11043 byte(s)
Added string comparison to all Value object classes (to compare the
string value of an object to another string).

Overloaded perl '.' operator to do dot product when the operands are
formulas returning vectors.  (Part of the auto-generation of
formulas).

A few improvements to real and complex class output results.

Made Union class slightly more robust and removed need for makeUnion
method other than in the Union itself.

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

aubreyja at gmail dot com
ViewVC Help
Powered by ViewVC 1.0.9