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

Diff of /trunk/pg/lib/Matrix.pm

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

Revision 1079 Revision 1280
4 4
5} 5}
6package Matrix; 6package Matrix;
7@Matrix::ISA = qw(MatrixReal1); 7@Matrix::ISA = qw(MatrixReal1);
8 8
9 9use Carp;
10 10
11$Matrix::DEFAULT_FORMAT = '% #-19.12E '; 11$Matrix::DEFAULT_FORMAT = '% #-19.12E ';
12# allows specification of the format 12# allows specification of the format
13sub _stringify 13sub _stringify
14{ 14{
172 } 172 }
173 } 173 }
174 } 174 }
175 return($matrix); 175 return($matrix);
176} 176}
177######################################################################
178# Modifications to MatrixReal.pm which allow use of complex entries
179######################################################################
180
181sub cp { # MEG makes new copies of complex number
182 my $z = shift;
183 return $z unless ref($z);
184 my $w = Complex1::cplx($z->Re,$z->Im);
185 return $w;
186}
187sub copy
188{
189 croak "Usage: \$matrix1->copy(\$matrix2);"
190 if (@_ != 2);
191
192 my($matrix1,$matrix2) = @_;
193 my($rows1,$cols1) = ($matrix1->[1],$matrix1->[2]);
194 my($rows2,$cols2) = ($matrix2->[1],$matrix2->[2]);
195 my($i,$j);
196
197 croak "MatrixReal1::copy(): matrix size mismatch"
198 unless (($rows1 == $rows2) && ($cols1 == $cols2));
199
200 for ( $i = 0; $i < $rows1; $i++ )
201 {
202 my $r1 = []; # New array ref
203 my $r2 = $matrix2->[0][$i];
204 #@$r1 = @$r2; # Copy whole array directly #MEG
205 # if the array contains complex objects new objects must be created.
206 foreach (@$r2) {
207 push(@$r1,cp($_) );
208 }
209 $matrix1->[0][$i] = $r1;
210 }
211 $matrix1->[3] = $matrix2->[3]; # sign or option
212 if (defined $matrix2->[4]) # is an LR decomposition matrix!
213 {
214 # $matrix1->[3] = $matrix2->[3]; # $sign
215 $matrix1->[4] = $matrix2->[4]; # $perm_row
216 $matrix1->[5] = $matrix2->[5]; # $perm_col
217 $matrix1->[6] = $matrix2->[6]; # $option
218 }
219}
220###################################################################
221
222# MEG added 6/25/03 to accomodate complex entries
223sub conj {
224 my $elem = shift;
225 $elem = (ref($elem)) ? ($elem->conjugate) : $elem;
226 $elem;
227}
228sub transpose
229{
230 croak "Usage: \$matrix1->transpose(\$matrix2);"
231 if (@_ != 2);
232
233 my($matrix1,$matrix2) = @_;
234 my($rows1,$cols1) = ($matrix1->[1],$matrix1->[2]);
235 my($rows2,$cols2) = ($matrix2->[1],$matrix2->[2]);
236
237 croak "MatrixReal1::transpose(): matrix size mismatch"
238 unless (($rows1 == $cols2) && ($cols1 == $rows2));
239
240 $matrix1->_undo_LR();
241
242 if ($rows1 == $cols1)
243 {
244 # more complicated to make in-place possible!
245 # # conj added by MEG
246 for (my $i = 0; $i < $rows1; $i++)
247 {
248 for (my $j = ($i + 1); $j < $cols1; $j++)
249 {
250 my $swap = conj($matrix2->[0][$i][$j]);
251 $matrix1->[0][$i][$j] = conj($matrix2->[0][$j][$i]);
252 $matrix1->[0][$j][$i] = $swap;
253 }
254 $matrix1->[0][$i][$i] = conj($matrix2->[0][$i][$i]);
255 }
256
257 }
258 else # ($rows1 != $cols1)
259 {
260 for (my $i = 0; $i < $rows1; $i++)
261 {
262 for (my $j = 0; $j < $cols1; $j++)
263 {
264 $matrix1->[0][$i][$j] = conj($matrix2->[0][$j][$i]);
265 }
266 }
267 }
268 $matrix1;
269}
270
177 271
1781; 2721;

Legend:
Removed from v.1079  
changed lines
  Added in v.1280

aubreyja at gmail dot com
ViewVC Help
Powered by ViewVC 1.0.9