Parent Directory
|
Revision Log
Revision 2579 - (view) (download) (as text)
| 1 : | sh002i | 2558 | ########################################################################### |
| 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' => \&compare, | ||
| 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->matrixMatrix(@{$M}) if (ref($M->[0]) eq 'ARRAY'); | ||
| 41 : | return $self->numberMatrix(@{$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 $d = scalar(@{$_[0]}); 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 : | dpvc | 2579 | $x = Value::Real->make($x) if $$Value::context->flag('useFuzzyReals') && !Value::isFormula($x); |
| 77 : | sh002i | 2558 | 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 : | dpvc | 2579 | my $open = $$Value::context->lists->get('Matrix')->{open}; |
| 341 : | my $close = $$Value::context->lists->get('Matrix')->{close}; | ||
| 342 : | sh002i | 2558 | return $open.join(',',@{$self->data}).$close |
| 343 : | if (Value::class($self->data->[0]) ne 'Matrix'); | ||
| 344 : | return $open.join(",\n ",@{$self->data}).$close; | ||
| 345 : | } | ||
| 346 : | |||
| 347 : | sub string { | ||
| 348 : | my $self = shift; my $equation = shift; | ||
| 349 : | dpvc | 2579 | my $open = shift || $$Value::context->lists->get('Matrix')->{open}; |
| 350 : | my $close = shift || $$Value::context->lists->get('Matrix')->{close}; | ||
| 351 : | sh002i | 2558 | my @coords = (); |
| 352 : | foreach my $x (@{$self->data}) { | ||
| 353 : | if (Value::isValue($x)) {push(@coords,$x->string($equation,$open,$close))} | ||
| 354 : | else {push(@coords,$x)} | ||
| 355 : | } | ||
| 356 : | return $open.join(',',@coords).$close; | ||
| 357 : | } | ||
| 358 : | |||
| 359 : | # | ||
| 360 : | # Use \matrix to lay out matrices | ||
| 361 : | # | ||
| 362 : | sub TeX { | ||
| 363 : | my $self = shift; my $equation = shift; | ||
| 364 : | dpvc | 2579 | my $open = shift || $$Value::context->lists->get('Matrix')->{open}; |
| 365 : | my $close = shift || $$Value::context->lists->get('Matrix')->{close}; | ||
| 366 : | sh002i | 2558 | $open = '\{' if $open eq '{'; $close = '\}' if $close eq '}'; |
| 367 : | my $TeX = ''; my @entries = (); | ||
| 368 : | dpvc | 2579 | if ($self->isRow) { |
| 369 : | foreach my $x (@{$self->data}) { | ||
| 370 : | sh002i | 2558 | push(@entries,(Value::isValue($x))? $x->TeX($equation,$open,$close): $x); |
| 371 : | } | ||
| 372 : | dpvc | 2579 | $TeX .= join(' &',@entries) . "\n"; |
| 373 : | } else { | ||
| 374 : | foreach my $row (@{$self->data}) { | ||
| 375 : | foreach my $x (@{$row->data}) { | ||
| 376 : | push(@entries,(Value::isValue($x))? $x->TeX($equation,$open,$close): $x); | ||
| 377 : | } | ||
| 378 : | $TeX .= join(' &',@entries) . '\cr'."\n"; @entries = (); | ||
| 379 : | } | ||
| 380 : | sh002i | 2558 | } |
| 381 : | return '\left'.$open.'\matrix{'."\n".$TeX.'}\right'.$close; | ||
| 382 : | } | ||
| 383 : | |||
| 384 : | ########################################################################### | ||
| 385 : | |||
| 386 : | 1; | ||
| 387 : |
| aubreyja at gmail dot com | ViewVC Help |
| Powered by ViewVC 1.0.9 |