source: src/main/java/agents/anac/y2015/agentBuyogV2/flanagan/complex/ComplexMatrix.java

Last change on this file was 127, checked in by Wouter Pasman, 6 years ago

#41 ROLL BACK of rev.126 . So this version is equal to rev. 125

File size: 53.4 KB
Line 
1/*
2* Class ComplexMatrix
3*
4* Defines a complex matrix and includes the methods
5* needed for standard matrix manipulations, e.g. multiplation,
6* and related procedures, e.g. solution of complex linear
7* simultaneous equations
8*
9* See class PhasorMatrix for phasor matrix manipulation
10* See class Complex for standard complex arithmetic
11*
12* WRITTEN BY: Dr Michael Thomas Flanagan
13*
14* DATE: June 2002
15* UPDATES: 16 February 2006, 7 March 2006, 31 March 2006, 22 April 2006, 1 July 2007, 9 April 2008
16* 18 October 2008, 16 June 2009, November 2009, 12 January 2010, 14 November 2010
17*
18*
19* DOCUMENTATION:
20* See Michael Thomas Flanagan's Java library on-line web pages:
21* http://www.ee.ucl.ac.uk/~mflanaga/java/ComplexMatrix.html
22* http://www.ee.ucl.ac.uk/~mflanaga/java/
23*
24* Copyright (c) 2002 - 2010 Michael Thomas Flanagan
25
26*
27* PERMISSION TO COPY:
28* Permission to use, copy and modify this software and its documentation for
29* NON-COMMERCIAL purposes is granted, without fee, provided that an acknowledgement
30* to the author, Michael Thomas Flanagan at www.ee.ucl.ac.uk/~mflanaga, appears in all copies.
31*
32* Dr Michael Thomas Flanagan makes no representations about the suitability
33* or fitness of the software for any or for a particular purpose.
34* Michael Thomas Flanagan shall not be liable for any damages suffered
35* as a result of using, modifying or distributing this software or its derivatives.
36*
37***************************************************************************************/
38
39package agents.anac.y2015.agentBuyogV2.flanagan.complex;
40
41import agents.anac.y2015.agentBuyogV2.flanagan.math.*;
42
43public class ComplexMatrix{
44
45 private int nrow = 0; // number of rows
46 private int ncol = 0; // number of columns
47 private Complex matrix[][] = null; // 2-D Complex Matrix
48 private int index[] = null; // row permutation index
49 private double dswap = 1.0D; // row swap index
50 private static final double TINY = 1.0e-30;
51
52/*********************************************************/
53
54 // CONSTRUCTORS
55 // Construct a nrow x ncol matrix of complex variables all equal to zero
56 public ComplexMatrix(int nrow, int ncol){
57 this.nrow = nrow;
58 this.ncol = ncol;
59 this.matrix = Complex.twoDarray(nrow, ncol);
60 this.index = new int[nrow];
61 for(int i=0;i<nrow;i++)this.index[i]=i;
62 this.dswap=1.0;
63 }
64
65 // Construct a nrow x ncol matrix of complex variables all equal to the complex number const
66 public ComplexMatrix(int nrow, int ncol, Complex constant){
67 this.nrow = nrow;
68 this.ncol = ncol;
69 this.matrix = Complex.twoDarray(nrow, ncol, constant);
70 this.index = new int[nrow];
71 for(int i=0;i<nrow;i++)this.index[i]=i;
72 this.dswap=1.0;
73 }
74
75 // Construct matrix with a copy of an existing nrow x ncol 2-D array of complex variables
76 public ComplexMatrix(Complex[][] twoD){
77 this.nrow = twoD.length;
78 this.ncol = twoD[0].length;
79 this.matrix = Complex.twoDarray(nrow, ncol);
80 for(int i=0; i<nrow; i++){
81 if(twoD[i].length!=ncol)throw new IllegalArgumentException("All rows must have the same length");
82 for(int j=0; j<ncol; j++){
83 this.matrix[i][j]=Complex.copy(twoD[i][j]);
84 }
85 }
86 this.index = new int[nrow];
87 for(int i=0;i<nrow;i++)this.index[i]=i;
88 this.dswap=1.0;
89 }
90
91 // Construct matrix with a copy of an existing nrow x ncol 2-D array of double variables
92 public ComplexMatrix(double[][] twoD){
93 this.nrow = twoD.length;
94 this.ncol = twoD[0].length;
95 for(int i=0; i<nrow; i++){
96 if(twoD[i].length!=ncol)throw new IllegalArgumentException("All rows must have the same length");
97 }
98 this.matrix = Complex.twoDarray(nrow, ncol);
99 for(int i=0; i<nrow; i++){
100 for(int j=0; j<ncol; j++){
101 this.matrix[i][j] = new Complex(twoD[i][j], 0.0);
102 }
103 }
104 this.index = new int[nrow];
105 for(int i=0;i<nrow;i++)this.index[i]=i;
106 this.dswap=1.0;
107 }
108
109 // Construct matrix with a copy of the complex matrix and permutation index of an existing ComplexMatrix bb.
110 public ComplexMatrix(ComplexMatrix bb){
111 this.nrow = bb.nrow;
112 this.ncol = bb.ncol;
113 this.matrix = (bb.copy()).matrix;
114 this.index = bb.index;
115 this.dswap = bb.dswap;
116 }
117
118 // Construct matrix with a copy of the 2D matrix and permutation index of an existing Matrix bb.
119 public ComplexMatrix(Matrix bb){
120 this.nrow = bb.getNrow();
121 this.ncol = bb.getNcol();
122 double[][] array = bb.getArrayCopy();
123 this.matrix = Complex.twoDarray(nrow, ncol);
124 for(int i=0; i<nrow; i++){
125 for(int j=0; j<ncol; j++){
126 this.matrix[i][j] = new Complex(array[i][j], 0.0);
127 }
128 }
129 this.index = bb.getIndexCopy();
130 this.dswap = bb.getSwap();
131 }
132
133
134 // SET VALUES
135 // Set the matrix with a copy of an existing nrow x ncol 2-D matrix of complex variables
136 public void setTwoDarray(Complex[][] aarray){
137 if(this.nrow != aarray.length)throw new IllegalArgumentException("row length of this ComplexMatrix differs from that of the 2D array argument");
138 if(this.ncol != aarray[0].length)throw new IllegalArgumentException("column length of this ComplexMatrix differs from that of the 2D array argument");
139 for(int i=0; i<nrow; i++){
140 if(aarray[i].length!=ncol)throw new IllegalArgumentException("All rows must have the same length");
141 for(int j=0; j<ncol; j++){
142 this.matrix[i][j]=Complex.copy(aarray[i][j]);
143 }
144 }
145 }
146
147 // Set the matrix with a copy of an existing nrow x ncol 2-D matrix of double variables
148 public void setTwoDarray(double[][] aarray){
149 if(this.nrow != aarray.length)throw new IllegalArgumentException("row length of this ComplexMatrix differs from that of the 2D array argument");
150 if(this.ncol != aarray[0].length)throw new IllegalArgumentException("column length of this ComplexMatrix differs from that of the 2D array argument");
151 for(int i=0; i<nrow; i++){
152 if(aarray[i].length!=ncol)throw new IllegalArgumentException("All rows must have the same length");
153 for(int j=0; j<ncol; j++){
154 this.matrix[i][j]=new Complex(aarray[i][j]);
155 }
156 }
157 }
158
159 // Set an individual array element
160 // i = row index
161 // j = column index
162 // aa = value of the element
163 public void setElement(int i, int j, Complex aa){
164 this.matrix[i][j]=Complex.copy(aa);
165 }
166
167 // Set an individual array element
168 // i = row index
169 // j = column index
170 // aa = real part of the element
171 // bb = imag part of the element
172 public void setElement(int i, int j, double aa, double bb){
173 this.matrix[i][j].reset(aa, bb);
174 }
175
176 // Set a sub-matrix starting with row index i, column index j
177 // and ending with row index k, column index l
178 // See SetSubMatrix below - this method has ben retained for compatibilityb purposes
179 public void setSubMatrix(int i, int j, int k, int l, Complex[][] subMatrix){
180 this.setSubMatrix(i, j, subMatrix);
181 }
182
183 // Set a sub-matrix starting with row index i, column index j
184 public void setSubMatrix(int i, int j, Complex[][] subMatrix){
185 int k = subMatrix.length;
186 int l = subMatrix[0].length;
187 if(i+k-1>=this.nrow)throw new IllegalArgumentException("Sub-matrix position is outside the row bounds of this Matrix");
188 if(j+l-1>=this.ncol)throw new IllegalArgumentException("Sub-matrix position is outside the column bounds of this Matrix");
189
190 int m = 0;
191 int n = 0;
192 for(int p=0; p<k; p++){
193 n = 0;
194 for(int q=0; q<l; q++){
195 this.matrix[i+p][j+q] = Complex.copy(subMatrix[m][n]);
196 n++;
197 }
198 m++;
199 }
200 }
201
202
203 // Set a sub-matrix
204 // row = array of row indices
205 // col = array of column indices
206 public void setSubMatrix(int[] row, int[] col, Complex[][] subMatrix){
207 int n=row.length;
208 int m=col.length;
209 for(int p=0; p<n; p++){
210 for(int q=0; q<m; q++){
211 this.matrix[row[p]][col[q]] = Complex.copy(subMatrix[p][q]);
212 }
213 }
214 }
215
216
217 // SPECIAL MATRICES
218 // Construct a complex identity matrix
219 public static ComplexMatrix identityMatrix(int nrow){
220 ComplexMatrix u = new ComplexMatrix(nrow, nrow);
221 for(int i=0; i<nrow; i++){
222 u.matrix[i][i]=Complex.plusOne();
223 }
224 return u;
225 }
226
227 // Construct a complex scalar matrix
228 public static ComplexMatrix scalarMatrix(int nrow, Complex diagconst){
229 ComplexMatrix u = new ComplexMatrix(nrow, nrow);
230 Complex[][] uarray = u.getArrayReference();
231 for(int i=0; i<nrow; i++){
232 for(int j=i; j<nrow; j++){
233 if(i==j){
234 uarray[i][j]=Complex.copy(diagconst);
235 }
236 }
237 }
238 return u;
239 }
240
241 // Construct a complex diagonal matrix
242 public static ComplexMatrix diagonalMatrix(int nrow, Complex[] diag){
243 if(diag.length!=nrow)throw new IllegalArgumentException("matrix dimension differs from diagonal array length");
244 ComplexMatrix u = new ComplexMatrix(nrow, nrow);
245 Complex[][] uarray = u.getArrayReference();
246 for(int i=0; i<nrow; i++){
247 for(int j=i; j<nrow; j++){
248 if(i==j){
249 uarray[i][j]=Complex.copy(diag[i]);
250 }
251 }
252 }
253 return u;
254 }
255
256 // COLUMN MATRICES
257 // Converts a 1-D array of Complex to a column matrix
258 public static ComplexMatrix columnMatrix(Complex[] darray){
259 int nr = darray.length;
260 ComplexMatrix pp = new ComplexMatrix(nr, 1);
261 for(int i=0; i<nr; i++)pp.matrix[i][0] = darray[i];
262 return pp;
263 }
264
265 // ROW MATRICES
266 // Converts a 1-D array of Complex to a row matrix
267 public static ComplexMatrix rowMatrix(Complex[] darray){
268 int nc = darray.length;
269 ComplexMatrix pp = new ComplexMatrix(1, nc);
270 for(int i=0; i<nc; i++)pp.matrix[0][i] = darray[i];
271 return pp;
272 }
273
274 // CONVERSIONS
275 // Converts a 1-D array of Complex to a complex column matrix
276 public static ComplexMatrix toComplexColumnMatrix(Complex[] carray){
277 int nr = carray.length;
278 ComplexMatrix cc = new ComplexMatrix(nr, 1);
279 for(int i=0; i<nr; i++)cc.matrix[i][0] = carray[i].copy();
280 return cc;
281 }
282
283 // Converts a 1-D array of doubles to a complex coumn matrix
284 public static ComplexMatrix toComplexColumnMatrix(double[] darray){
285 int nr = darray.length;
286 ComplexMatrix cc = new ComplexMatrix(nr, 1);
287 for(int i=0; i<nr; i++)cc.matrix[i][0].reset(darray[i], 0.0D);
288 return cc;
289 }
290
291 // Converts a 1-D array of Complex to a complex row matrix
292 public static ComplexMatrix toComplexRowMatrix(Complex[] carray){
293 int nc = carray.length;
294 ComplexMatrix cc = new ComplexMatrix(1, nc);
295 for(int i=0; i<nc; i++)cc.matrix[0][i] = carray[i].copy();
296 return cc;
297 }
298
299 // Converts a 1-D array of doubles to a complex row matrix
300 public static ComplexMatrix toComplexRowMatrix(double[] darray){
301 int nc = darray.length;
302 ComplexMatrix cc = new ComplexMatrix(1, nc);
303 for(int i=0; i<nc; i++)cc.matrix[0][i].reset(darray[i], 0.0D);
304 return cc;
305 }
306
307 // Converts a matrix of doubles (Matrix) to a complex matrix (ComplexMatix)
308 public static ComplexMatrix toComplexMatrix(Matrix marray){
309 int nr = marray.getNrow();
310 int nc = marray.getNcol();
311
312 ComplexMatrix pp = new ComplexMatrix(nr, nc);
313 for(int i=0; i<nr; i++){
314 for(int j=0; j<nc; j++){
315 pp.matrix[i][j].reset(marray.getElementCopy(i, j), 0.0D);
316 }
317 }
318 return pp;
319 }
320
321 // Converts a 2D array of doubles to a complex matrix (ComplexMatix)
322 public static ComplexMatrix toComplexMatrix(double[][] darray){
323 int nr = darray.length;
324 int nc = darray[0].length;
325 for(int i=1; i<nr; i++){
326 if(darray[i].length!=nc)throw new IllegalArgumentException("All rows must have the same length");
327 }
328 ComplexMatrix pp = new ComplexMatrix(nr, nc);
329 for(int i=0; i<pp.nrow; i++){
330 for(int j=0; j<pp.ncol; j++){
331 pp.matrix[i][j].reset(darray[i][j], 0.0D);
332 }
333 }
334 return pp;
335 }
336
337 // GET VALUES
338 // Return the number of rows
339 public int getNrow(){
340 return this.nrow;
341 }
342
343 // Return the number of columns
344 public int getNcol(){
345 return this.ncol;
346 }
347
348 // Return a reference to the internal 2-D array
349 public Complex[][] getArrayReference(){
350 return this.matrix;
351 }
352
353 // Return a reference to the internal 2-D array
354 public Complex[][] getArray(){
355 return this.matrix;
356 }
357
358 // Return a reference to the internal 2-D array
359 // included for backward compatibility with earlier incorrect documentation
360 public Complex[][] getArrayPointer(){
361 return this.matrix;
362 }
363
364 // Return a copy of the internal 2-D array
365 public Complex[][] getArrayCopy(){
366 Complex[][] c = new Complex[this.nrow][this.ncol];
367 for(int i=0; i<nrow; i++){
368 for(int j=0; j<ncol; j++){
369 c[i][j]=Complex.copy(matrix[i][j]);
370 }
371 }
372 return c;
373 }
374
375 // Return a single element of the internal 2-D array
376 public Complex getElementReference(int i, int j){
377 return this.matrix[i][j];
378 }
379
380 // Return a reference to a single element of the internal 2-D array
381 // included for backward compatibility with earlier incorrect documentation
382 public Complex getElementPointer(int i, int j){
383 return this.matrix[i][j];
384 }
385
386 // Return a copy of a single element of the internal 2-D array
387 public Complex getElementCopy(int i, int j){
388 return Complex.copy(this.matrix[i][j]);
389 }
390
391 // Return a sub-matrix starting with row index i, column index j
392 // and ending with row index k, column index l
393 public ComplexMatrix getSubMatrix(int i, int j, int k, int l){
394 if(i>k)throw new IllegalArgumentException("row indices inverted");
395 if(j>l)throw new IllegalArgumentException("column indices inverted");
396 if(k>=this.nrow)throw new IllegalArgumentException("Sub-matrix position is outside the row bounds of this Matrix" );
397 if(l>=this.ncol)throw new IllegalArgumentException("Sub-matrix position is outside the column bounds of this Matrix" + i + " " +l);
398
399 int n=k-i+1, m=l-j+1;
400 ComplexMatrix subMatrix = new ComplexMatrix(n, m);
401 Complex[][] sarray = subMatrix.getArrayReference();
402 for(int p=0; p<n; p++){
403 for(int q=0; q<m; q++){
404 sarray[p][q]=Complex.copy(this.matrix[i+p][j+q]);
405 }
406 }
407 return subMatrix;
408 }
409
410 // Return a sub-matrix
411 // row = array of row indices
412 // col = array of column indices
413 public ComplexMatrix getSubMatrix(int[] row, int[] col){
414 int n = row.length;
415 int m = col.length;
416 ComplexMatrix subMatrix = new ComplexMatrix(n, m);
417 Complex[][] sarray = subMatrix.getArrayReference();
418 for(int i=0; i<n; i++){
419 for(int j=0; j<m; j++){
420 sarray[i][j]=Complex.copy(this.matrix[row[i]][col[j]]);
421 }
422 }
423 return subMatrix;
424 }
425
426 // Return a reference to the permutation index array
427 public int[] getIndexReference(){
428 return this.index;
429 }
430
431 // Return a reference to the permutation index array
432 public int[] getIndexPointer(){
433 return this.index;
434 }
435
436 // Return a copy of the permutation index array
437 public int[] getIndexCopy(){
438 int[] indcopy = new int[this.nrow];
439 for(int i=0; i<this.nrow; i++){
440 indcopy[i]=this.index[i];
441 }
442 return indcopy;
443 }
444
445 // Return the row swap index
446 public double getSwap(){
447 return this.dswap;
448 }
449
450 // COPY
451 // Copy a ComplexMatrix [static method]
452 public static ComplexMatrix copy(ComplexMatrix a){
453 if(a==null){
454 return null;
455 }
456 else{
457 int nr = a.getNrow();
458 int nc = a.getNcol();
459 Complex[][] aarray = a.getArrayReference();
460 ComplexMatrix b = new ComplexMatrix(nr,nc);
461 b.nrow = nr;
462 b.ncol = nc;
463 Complex[][] barray = b.getArrayReference();
464 for(int i=0; i<nr; i++){
465 for(int j=0; j<nc; j++){
466 barray[i][j]=Complex.copy(aarray[i][j]);
467 }
468 }
469 for(int i=0; i<nr; i++)b.index[i] = a.index[i];
470 return b;
471 }
472 }
473
474 // Copy a ComplexMatrix [instance method]
475 public ComplexMatrix copy(){
476 if(this==null){
477 return null;
478 }
479 else{
480 int nr = this.nrow;
481 int nc = this.ncol;
482 ComplexMatrix b = new ComplexMatrix(nr,nc);
483 Complex[][] barray = b.getArrayReference();
484 b.nrow = nr;
485 b.ncol = nc;
486 for(int i=0; i<nr; i++){
487 for(int j=0; j<nc; j++){
488 barray[i][j]=Complex.copy(this.matrix[i][j]);
489 }
490 }
491 for(int i=0; i<nr; i++)b.index[i] = this.index[i];
492 return b;
493 }
494 }
495
496 // Clone a ComplexMatrix
497 public Object clone(){
498 if(this==null){
499 return null;
500 }
501 else{
502 int nr = this.nrow;
503 int nc = this.ncol;
504 ComplexMatrix b = new ComplexMatrix(nr,nc);
505 Complex[][] barray = b.getArrayReference();
506 b.nrow = nr;
507 b.ncol = nc;
508 for(int i=0; i<nr; i++){
509 for(int j=0; j<nc; j++){
510 barray[i][j]=Complex.copy(this.matrix[i][j]);
511 }
512 }
513 for(int i=0; i<nr; i++)b.index[i] = this.index[i];
514 return (Object) b;
515 }
516 }
517
518 // ADDITION
519 // Add this matrix to matrix B. This matrix remains unaltered [instance method]
520 public ComplexMatrix plus(ComplexMatrix bmat){
521 if((this.nrow!=bmat.nrow)||(this.ncol!=bmat.ncol)){
522 throw new IllegalArgumentException("Array dimensions do not agree");
523 }
524 int nr=bmat.nrow;
525 int nc=bmat.ncol;
526 ComplexMatrix cmat = new ComplexMatrix(nr,nc);
527 Complex[][] carray = cmat.getArrayReference();
528 for(int i=0; i<nr; i++){
529 for(int j=0; j<nc; j++){
530 carray[i][j]=this.matrix[i][j].plus(bmat.matrix[i][j]);
531 }
532 }
533 return cmat;
534 }
535
536 // Add this matrix to a Comlex 2-D array. [instance method]
537 public ComplexMatrix plus(Complex[][] bmat){
538 int nr=bmat.length;
539 int nc=bmat[0].length;
540 if((this.nrow!=nr)||(this.ncol!=nc)){
541 throw new IllegalArgumentException("Array dimensions do not agree");
542 }
543 ComplexMatrix cmat = new ComplexMatrix(nr,nc);
544 Complex[][] carray = cmat.getArrayReference();
545 for(int i=0; i<nr; i++){
546 for(int j=0; j<nc; j++){
547 carray[i][j]=this.matrix[i][j].plus(bmat[i][j]);
548 }
549 }
550 return cmat;
551 }
552
553 // Add this matrix to a real matrix B. [instance method]
554 public ComplexMatrix plus(Matrix bmat){
555 int nr=bmat.getNrow();
556 int nc=bmat.getNcol();
557 if((this.nrow!=nr)||(this.ncol!=nc)){
558 throw new IllegalArgumentException("Array dimensions do not agree");
559 }
560
561 ComplexMatrix cmat = new ComplexMatrix(nr,nc);
562 Complex[][] carray = cmat.getArrayReference();
563 for(int i=0; i<nr; i++){
564 for(int j=0; j<nc; j++){
565 carray[i][j]=this.matrix[i][j].plus(bmat.getElement(i,j));
566 }
567 }
568 return cmat;
569 }
570
571 // Add this matrix to a real 2-D array. [instance method]
572 public ComplexMatrix plus(double[][] bmat){
573 int nr=bmat.length;
574 int nc=bmat[0].length;
575 if((this.nrow!=nr)||(this.ncol!=nc)){
576 throw new IllegalArgumentException("Array dimensions do not agree");
577 }
578
579 ComplexMatrix cmat = new ComplexMatrix(nr,nc);
580 Complex[][] carray = cmat.getArrayReference();
581 for(int i=0; i<nr; i++){
582 for(int j=0; j<nc; j++){
583 carray[i][j]=this.matrix[i][j].plus(bmat[i][j]);
584 }
585 }
586 return cmat;
587 }
588
589 // Add matrices A and B [static method]
590 public static ComplexMatrix plus(ComplexMatrix amat, ComplexMatrix bmat){
591 if((amat.nrow!=bmat.nrow)||(amat.ncol!=bmat.ncol)){
592 throw new IllegalArgumentException("Array dimensions do not agree");
593 }
594 int nr=amat.nrow;
595 int nc=amat.ncol;
596 ComplexMatrix cmat = new ComplexMatrix(nr,nc);
597 Complex[][] carray = cmat.getArrayReference();
598 for(int i=0; i<nr; i++){
599 for(int j=0; j<nc; j++){
600 carray[i][j]=amat.matrix[i][j].plus(bmat.matrix[i][j]);
601 }
602 }
603 return cmat;
604 }
605
606 // Add matrix B to this matrix [equivalence of +=]
607 public void plusEquals(ComplexMatrix bmat){
608 if((this.nrow!=bmat.nrow)||(this.ncol!=bmat.ncol)){
609 throw new IllegalArgumentException("Array dimensions do not agree");
610 }
611 int nr=bmat.nrow;
612 int nc=bmat.ncol;
613
614 for(int i=0; i<nr; i++){
615 for(int j=0; j<nc; j++){
616 this.matrix[i][j].plusEquals(bmat.matrix[i][j]);
617 }
618 }
619 }
620
621 // SUBTRACTION
622 // Subtract matrix B from this matrix. This matrix remains unaltered [instance method]
623 public ComplexMatrix minus(ComplexMatrix bmat){
624 if((this.nrow!=bmat.nrow)||(this.ncol!=bmat.ncol)){
625 throw new IllegalArgumentException("Array dimensions do not agree");
626 }
627 int nr=this.nrow;
628 int nc=this.ncol;
629 ComplexMatrix cmat = new ComplexMatrix(nr,nc);
630 Complex[][] carray = cmat.getArrayReference();
631 for(int i=0; i<nr; i++){
632 for(int j=0; j<nc; j++){
633 carray[i][j]=this.matrix[i][j].minus(bmat.matrix[i][j]);
634 }
635 }
636 return cmat;
637 }
638
639 // Subtract Comlex 2-D array from this matrix. [instance method]
640 public ComplexMatrix minus(Complex[][] bmat){
641 int nr=bmat.length;
642 int nc=bmat[0].length;
643 if((this.nrow!=nr)||(this.ncol!=nc)){
644 throw new IllegalArgumentException("Array dimensions do not agree");
645 }
646 ComplexMatrix cmat = new ComplexMatrix(nr,nc);
647 Complex[][] carray = cmat.getArrayReference();
648 for(int i=0; i<nr; i++){
649 for(int j=0; j<nc; j++){
650 carray[i][j]=this.matrix[i][j].minus(bmat[i][j]);
651 }
652 }
653 return cmat;
654 }
655
656 // Subtract a real matrix from a real matrix B. [instance method]
657 public ComplexMatrix minus(Matrix bmat){
658 int nr=bmat.getNrow();
659 int nc=bmat.getNcol();
660 if((this.nrow!=nr)||(this.ncol!=nc)){
661 throw new IllegalArgumentException("Array dimensions do not agree");
662 }
663
664 ComplexMatrix cmat = new ComplexMatrix(nr,nc);
665 Complex[][] carray = cmat.getArrayReference();
666 for(int i=0; i<nr; i++){
667 for(int j=0; j<nc; j++){
668 carray[i][j]=this.matrix[i][j].minus(bmat.getElement(i,j));
669 }
670 }
671 return cmat;
672 }
673
674 // Subtract a real 2-D array from this matrix. [instance method]
675 public ComplexMatrix minus(double[][] bmat){
676 int nr=bmat.length;
677 int nc=bmat[0].length;
678 if((this.nrow!=nr)||(this.ncol!=nc)){
679 throw new IllegalArgumentException("Array dimensions do not agree");
680 }
681
682 ComplexMatrix cmat = new ComplexMatrix(nr,nc);
683 Complex[][] carray = cmat.getArrayReference();
684 for(int i=0; i<nr; i++){
685 for(int j=0; j<nc; j++){
686 carray[i][j]=this.matrix[i][j].minus(bmat[i][j]);
687 }
688 }
689 return cmat;
690 }
691
692
693 // Subtract matrix B from matrix A [static method]
694 public static ComplexMatrix minus(ComplexMatrix amat, ComplexMatrix bmat){
695 if((amat.nrow!=bmat.nrow)||(amat.ncol!=bmat.ncol)){
696 throw new IllegalArgumentException("Array dimensions do not agree");
697 }
698 int nr=amat.nrow;
699 int nc=amat.ncol;
700 ComplexMatrix cmat = new ComplexMatrix(nr,nc);
701 Complex[][] carray = cmat.getArrayReference();
702 for(int i=0; i<nr; i++){
703 for(int j=0; j<nc; j++){
704 carray[i][j]=amat.matrix[i][j].minus(bmat.matrix[i][j]);
705 }
706 }
707 return cmat;
708 }
709
710 // Subtract matrix B from this matrix [equivlance of -=]
711 public void minusEquals(ComplexMatrix bmat){
712 if((this.nrow!=bmat.nrow)||(this.ncol!=bmat.ncol)){
713 throw new IllegalArgumentException("Array dimensions do not agree");
714 }
715 int nr=bmat.nrow;
716 int nc=bmat.ncol;
717
718 for(int i=0; i<nr; i++){
719 for(int j=0; j<nc; j++){
720 this.matrix[i][j].minusEquals(bmat.matrix[i][j]);
721 }
722 }
723 }
724
725 // MULTIPLICATION
726 // Multiply this complex matrix by a complex matrix. [instance method]
727 // This matrix remains unaltered.
728 public ComplexMatrix times(ComplexMatrix bmat){
729 if(this.ncol!=bmat.nrow)throw new IllegalArgumentException("Nonconformable matrices");
730
731 ComplexMatrix cmat = new ComplexMatrix(this.nrow, bmat.ncol);
732 Complex [][] carray = cmat.getArrayReference();
733 Complex sum = new Complex();
734
735 for(int i=0; i<this.nrow; i++){
736 for(int j=0; j<bmat.ncol; j++){
737 sum=Complex.zero();
738 for(int k=0; k<this.ncol; k++){
739 sum.plusEquals(this.matrix[i][k].times(bmat.matrix[k][j]));
740 }
741 carray[i][j]=Complex.copy(sum);
742 }
743 }
744 return cmat;
745 }
746
747 // Multiply this complex matrix by a complex 2-D array. [instance method]
748 public ComplexMatrix times(Complex[][] bmat){
749 int nr=bmat.length;
750 int nc=bmat[0].length;
751 if(this.ncol!=nr)throw new IllegalArgumentException("Nonconformable matrices");
752
753 ComplexMatrix cmat = new ComplexMatrix(this.nrow, nc);
754 Complex [][] carray = cmat.getArrayReference();
755 Complex sum = new Complex();
756
757 for(int i=0; i<this.nrow; i++){
758 for(int j=0; j<nc; j++){
759 sum=Complex.zero();
760 for(int k=0; k<this.ncol; k++){
761 sum.plusEquals(this.matrix[i][k].times(bmat[k][j]));
762 }
763 carray[i][j]=Complex.copy(sum);
764 }
765 }
766 return cmat;
767 }
768
769 // Multiply this complex matrix by a real matrix. [instance method]
770 // This matrix remains unaltered.
771 public ComplexMatrix times(Matrix bmat){
772 int nr=bmat.getNrow();
773 int nc=bmat.getNcol();
774
775 if(this.ncol!=nr)throw new IllegalArgumentException("Nonconformable matrices");
776
777 ComplexMatrix cmat = new ComplexMatrix(this.nrow, nc);
778 Complex [][] carray = cmat.getArrayReference();
779 Complex sum = new Complex();
780
781 for(int i=0; i<this.nrow; i++){
782 for(int j=0; j<nc; j++){
783 sum=Complex.zero();
784 for(int k=0; k<this.ncol; k++){
785 sum.plusEquals(this.matrix[i][k].times(bmat.getElement(k,j)));
786 }
787 carray[i][j]=Complex.copy(sum);
788 }
789 }
790 return cmat;
791 }
792
793 // Multiply this complex matrix by a real 2-D array. [instance method]
794 public ComplexMatrix times(double[][] bmat){
795 int nr=bmat.length;
796 int nc=bmat[0].length;
797 if(this.ncol!=nr)throw new IllegalArgumentException("Nonconformable matrices");
798
799 ComplexMatrix cmat = new ComplexMatrix(this.nrow, nc);
800 Complex [][] carray = cmat.getArrayReference();
801 Complex sum = new Complex();
802
803 for(int i=0; i<this.nrow; i++){
804 for(int j=0; j<nc; j++){
805 sum=Complex.zero();
806 for(int k=0; k<this.ncol; k++){
807 sum.plusEquals(this.matrix[i][k].times(bmat[k][j]));
808 }
809 carray[i][j]=Complex.copy(sum);
810 }
811 }
812 return cmat;
813 }
814
815 // Multiply this complex matrix by a complex constant [instance method]
816 // This matrix remains unaltered
817 public ComplexMatrix times(Complex constant){
818 ComplexMatrix cmat = new ComplexMatrix(this.nrow, this.ncol);
819 Complex [][] carray = cmat.getArrayReference();
820
821 for(int i=0; i<this.nrow; i++){
822 for(int j=0; j<this.ncol; j++){
823 carray[i][j] = this.matrix[i][j].times(constant);
824 }
825 }
826 return cmat;
827 }
828
829 // Multiply this complex matrix by a real (double) constant [instance method]
830 // This matrix remains unaltered.
831 public ComplexMatrix times(double constant){
832 ComplexMatrix cmat = new ComplexMatrix(this.nrow, this.ncol);
833 Complex [][] carray = cmat.getArrayReference();
834 Complex cconstant = new Complex(constant, 0.0);
835
836 for(int i=0; i<this.nrow; i++){
837 for(int j=0; j<this.ncol; j++){
838 carray[i][j] = this.matrix[i][j].times(cconstant);
839 }
840 }
841 return cmat;
842 }
843
844 // Multiply two complex matrices {static method]
845 public static ComplexMatrix times(ComplexMatrix amat, ComplexMatrix bmat){
846 if(amat.ncol!=bmat.nrow)throw new IllegalArgumentException("Nonconformable matrices");
847
848 ComplexMatrix cmat = new ComplexMatrix(amat.nrow, bmat.ncol);
849 Complex [][] carray = cmat.getArrayReference();
850 Complex sum = new Complex();
851
852 for(int i=0; i<amat.nrow; i++){
853 for(int j=0; j<bmat.ncol; j++){
854 sum=Complex.zero();
855 for(int k=0; k<amat.ncol; k++){
856 sum.plusEquals(amat.matrix[i][k].times(bmat.matrix[k][j]));
857 }
858 carray[i][j]=Complex.copy(sum);
859 }
860 }
861 return cmat;
862 }
863
864
865 // Multiply a complex matrix by a complex constant [static method]
866 public static ComplexMatrix times(ComplexMatrix amat, Complex constant){
867 ComplexMatrix cmat = new ComplexMatrix(amat.nrow, amat.ncol);
868 Complex [][] carray = cmat.getArrayReference();
869
870 for(int i=0; i<amat.nrow; i++){
871 for(int j=0; j<amat.ncol; j++){
872 carray[i][j] = amat.matrix[i][j].times(constant);
873 }
874 }
875 return cmat;
876 }
877
878 // Multiply a complex matrix by a real (double) constant [static method]
879 public static ComplexMatrix times(ComplexMatrix amat, double constant){
880 ComplexMatrix cmat = new ComplexMatrix(amat.nrow, amat.ncol);
881 Complex [][] carray = cmat.getArrayReference();
882 Complex cconstant = new Complex(constant, 0.0);
883
884 for(int i=0; i<amat.nrow; i++){
885 for(int j=0; j<amat.ncol; j++){
886 carray[i][j] = amat.matrix[i][j].times(cconstant);
887 }
888 }
889 return cmat;
890 }
891
892 // Multiply this matrix by a complex matrix [equivalence of *=]
893 public void timesEquals(ComplexMatrix bmat){
894 if(this.ncol!=bmat.nrow)throw new IllegalArgumentException("Nonconformable matrices");
895
896 ComplexMatrix cmat = new ComplexMatrix(this.nrow, bmat.ncol);
897 Complex [][] carray = cmat.getArrayReference();
898 Complex sum = new Complex();
899
900 for(int i=0; i<this.nrow; i++){
901 for(int j=0; j<bmat.ncol; j++){
902 sum=Complex.zero();
903 for(int k=0; k<this.ncol; k++){
904 sum.plusEquals(this.matrix[i][k].times(bmat.matrix[k][j]));
905 }
906 carray[i][j]=Complex.copy(sum);
907 }
908 }
909
910 this.nrow = cmat.nrow;
911 this.ncol = cmat.ncol;
912 for(int i=0; i<this.nrow; i++){
913 for(int j=0; j<this.ncol; j++){
914 this.matrix[i][j] = cmat.matrix[i][j];
915 }
916 }
917 }
918
919
920
921
922 // Multiply this matrix by a complex constant [equivalence of *=]
923 public void timesEquals(Complex constant){
924
925 for(int i=0; i<this.nrow; i++){
926 for(int j=0; j<this.ncol; j++){
927 this.matrix[i][j].timesEquals(constant);
928 }
929 }
930 }
931
932 // Multiply this matrix by a real (double) constant [equivalence of *=]
933 public void timesEquals(double constant){
934 Complex cconstant = new Complex(constant, 0.0);
935
936 for(int i=0; i<this.nrow; i++){
937 for(int j=0; j<this.ncol; j++){
938 this.matrix[i][j].timesEquals(cconstant);
939 }
940 }
941 }
942
943 // DIVISION
944 // Divide this ComplexMatrix by a ComplexMatrix - instance method.
945 public ComplexMatrix over(ComplexMatrix bmat){
946 if((this.nrow!=bmat.nrow)||(this.ncol!=bmat.ncol)){
947 throw new IllegalArgumentException("Array dimensions do not agree");
948 }
949 return this.times(bmat.inverse());
950 }
951
952 // Divide this matrix by a Complex 2-D array - instance method.
953 public ComplexMatrix over(Complex[][] bmat){
954 int nr=bmat.length;
955 int nc=bmat[0].length;
956 if((this.nrow!=nr)||(this.ncol!=nc)){
957 throw new IllegalArgumentException("Array dimensions do not agree");
958 }
959
960 ComplexMatrix cmat = new ComplexMatrix(bmat);
961 return this.times(cmat.inverse());
962 }
963
964 // Divide this ComplexMatrix by a Matrix - instance method.
965 public ComplexMatrix over(Matrix bmat){
966 ComplexMatrix pmat = ComplexMatrix.toComplexMatrix(bmat);
967 return this.over(pmat);
968 }
969
970 // Divide this ComplexMatrix by a 2D array of double - instance method.
971 public ComplexMatrix over(double[][] bmat){
972 ComplexMatrix pmat = ComplexMatrix.toComplexMatrix(bmat);
973 return this.over(pmat);
974 }
975
976 // Divide this ComplexMatrix by a ComplexMatrix - static method.
977 public ComplexMatrix over(ComplexMatrix amat, ComplexMatrix bmat){
978 if((amat.nrow!=bmat.nrow)||(amat.ncol!=bmat.ncol)){
979 throw new IllegalArgumentException("Array dimensions do not agree");
980 }
981 return amat.times(bmat.inverse());
982 }
983
984 // Divide this ComplexMatrix by a ComplexMatrix [equivalence of /=]
985 public void overEquals(ComplexMatrix bmat){
986 if((this.nrow!=bmat.nrow)||(this.ncol!=bmat.ncol)){
987 throw new IllegalArgumentException("Array dimensions do not agree");
988 }
989 ComplexMatrix cmat = new ComplexMatrix(bmat);
990 this.timesEquals(cmat.inverse());
991 }
992
993 // INVERSE
994 // Inverse of a square complex matrix [instance method]
995 public ComplexMatrix inverse(){
996 int n = this.nrow;
997 if(n!=this.ncol)throw new IllegalArgumentException("Matrix is not square");
998 ComplexMatrix invmat = new ComplexMatrix(n, n);
999
1000 if(n==1){
1001 Complex[][] hold = this.getArrayCopy();
1002 if(hold[0][0].isZero())throw new IllegalArgumentException("Matrix is singular");
1003 hold[0][0] = Complex.plusOne().over(hold[0][0]);
1004 invmat = new ComplexMatrix(hold);
1005 }
1006 else{
1007 if(n==2){
1008 Complex[][] hold = this.getArrayCopy();
1009 Complex det = (hold[0][0].times(hold[1][1])).minus(hold[0][1].times(hold[1][0]));
1010 if(det.isZero())throw new IllegalArgumentException("Matrix is singular");
1011
1012 Complex[][] hold2 = Complex.twoDarray(2,2);
1013 hold2[0][0] = hold[1][1].over(det);
1014 hold2[1][1] = hold[0][0].over(det);
1015 hold2[1][0] = hold[1][0].negate().over(det);
1016 hold2[0][1] = hold[0][1].negate().over(det);
1017 invmat = new ComplexMatrix(hold2);
1018 }
1019 else{
1020 Complex[] col = new Complex[n];
1021 Complex[] xvec = new Complex[n];
1022 Complex[][] invarray = invmat.getArrayReference();
1023 ComplexMatrix ludmat;
1024
1025 ludmat = this.luDecomp();
1026 for(int j=0; j<n; j++){
1027 for(int i=0; i<n; i++)col[i]=Complex.zero();
1028 col[j]=Complex.plusOne();
1029 xvec=ludmat.luBackSub(col);
1030 for(int i=0; i<n; i++)invarray[i][j]=Complex.copy(xvec[i]);
1031 }
1032 }
1033 }
1034 return invmat;
1035 }
1036
1037 // Inverse of a square complex matrix [static method]
1038 public static ComplexMatrix inverse(ComplexMatrix amat){
1039 int n = amat.nrow;
1040 if(n!=amat.ncol)throw new IllegalArgumentException("Matrix is not square");
1041
1042 ComplexMatrix invmat = new ComplexMatrix(n, n);
1043
1044 if(n==1){
1045 Complex[][] hold = amat.getArrayCopy();
1046 if(hold[0][0].isZero())throw new IllegalArgumentException("Matrix is singular");
1047 hold[0][0] = Complex.plusOne().over(hold[0][0]);
1048 invmat = new ComplexMatrix(hold);
1049 }
1050 else{
1051 if(n==2){
1052 Complex[][] hold = amat.getArrayCopy();
1053 Complex det = (hold[0][0].times(hold[1][1])).minus(hold[0][1].times(hold[1][0]));
1054 if(det.isZero())throw new IllegalArgumentException("Matrix is singular");
1055
1056 Complex[][] hold2 = Complex.twoDarray(2,2);
1057 hold2[0][0] = hold[1][1].over(det);
1058 hold2[1][1] = hold[0][0].over(det);
1059 hold2[1][0] = hold[1][0].negate().over(det);
1060 hold2[0][1] = hold[0][1].negate().over(det);
1061 invmat = new ComplexMatrix(hold2);
1062 }
1063 else{
1064 Complex[] col = new Complex[n];
1065 Complex[] xvec = new Complex[n];
1066 Complex[][] invarray = invmat.getArrayReference();
1067 ComplexMatrix ludmat;
1068
1069 ludmat = amat.luDecomp();
1070 for(int j=0; j<n; j++){
1071 for(int i=0; i<n; i++)col[i]=Complex.zero();
1072 col[j]=Complex.plusOne();
1073 xvec=ludmat.luBackSub(col);
1074 for(int i=0; i<n; i++)invarray[i][j]=Complex.copy(xvec[i]);
1075 }
1076 }
1077 }
1078 return invmat;
1079 }
1080
1081 // TRANSPOSE
1082 // Transpose of a complex matrix [instance method]
1083 public ComplexMatrix transpose(){
1084 ComplexMatrix tmat = new ComplexMatrix(this.ncol, this.nrow);
1085 Complex[][] tarray = tmat.getArrayReference();
1086 for(int i=0; i<this.ncol; i++){
1087 for(int j=0; j<this.nrow; j++){
1088 tarray[i][j]=Complex.copy(this.matrix[j][i]);
1089 }
1090 }
1091 return tmat;
1092 }
1093
1094 // Transpose of a complex matrix [static method]
1095 public static ComplexMatrix transpose(ComplexMatrix amat){
1096 ComplexMatrix tmat = new ComplexMatrix(amat.ncol, amat.nrow);
1097 Complex[][] tarray = tmat.getArrayReference();
1098 for(int i=0; i<amat.ncol; i++){
1099 for(int j=0; j<amat.nrow; j++){
1100 tarray[i][j]=Complex.copy(amat.matrix[j][i]);
1101 }
1102 }
1103 return tmat;
1104 }
1105
1106 // COMPLEX CONJUGATE
1107 //Complex Conjugate of a complex matrix [instance method]
1108 public ComplexMatrix conjugate(){
1109 ComplexMatrix conj = ComplexMatrix.copy(this);
1110 for(int i=0; i<this.nrow; i++){
1111 for(int j=0; j<this.ncol; j++){
1112 conj.matrix[i][j]=this.matrix[i][j].conjugate();
1113 }
1114 }
1115 return conj;
1116 }
1117
1118 //Complex Conjugate of a complex matrix [static method]
1119 public static ComplexMatrix conjugate(ComplexMatrix amat){
1120 ComplexMatrix conj = ComplexMatrix.copy(amat);
1121 for(int i=0; i<amat.nrow; i++){
1122 for(int j=0; j<amat.ncol; j++){
1123 conj.matrix[i][j]=amat.matrix[i][j].conjugate();
1124 }
1125 }
1126 return conj;
1127 }
1128
1129 // ADJOIN
1130 // Adjoin of a complex matrix [instance method]
1131 public ComplexMatrix adjoin(){
1132 ComplexMatrix adj = ComplexMatrix.copy(this);
1133 adj=adj.transpose();
1134 adj=adj.conjugate();
1135 return adj;
1136 }
1137
1138 // Adjoin of a complex matrix [static method]
1139 public ComplexMatrix adjoin(ComplexMatrix amat){
1140 ComplexMatrix adj = ComplexMatrix.copy(amat);
1141 adj=adj.transpose();
1142 adj=adj.conjugate();
1143 return adj;
1144 }
1145
1146 // OPPOSITE
1147 // Opposite of a complex matrix [instance method]
1148 public ComplexMatrix opposite(){
1149 ComplexMatrix opp = ComplexMatrix.copy(this);
1150 for(int i=0; i<this.nrow; i++){
1151 for(int j=0; j<this.ncol; j++){
1152 opp.matrix[i][j]=this.matrix[i][j].times(Complex.minusOne());
1153 }
1154 }
1155 return opp;
1156 }
1157
1158 // Opposite of a complex matrix [static method]
1159 public static ComplexMatrix opposite(ComplexMatrix amat){
1160 ComplexMatrix opp = ComplexMatrix.copy(amat);
1161 for(int i=0; i<amat.nrow; i++){
1162 for(int j=0; j<amat.ncol; j++){
1163 opp.matrix[i][j]=amat.matrix[i][j].times(Complex.minusOne());
1164 }
1165 }
1166 return opp;
1167 }
1168
1169 // TRACE
1170 // Trace of a complex matrix [instance method]
1171 public Complex trace(){
1172 Complex trac = new Complex(0.0, 0.0);
1173 for(int i=0; i<Math.min(this.ncol,this.ncol); i++){
1174 trac.plusEquals(this.matrix[i][i]);
1175 }
1176 return trac;
1177 }
1178
1179 // Trace of a complex matrix [static method]
1180 public static Complex trace(ComplexMatrix amat){
1181 Complex trac = new Complex(0.0, 0.0);
1182 for(int i=0; i<Math.min(amat.ncol,amat.ncol); i++){
1183 trac.plusEquals(amat.matrix[i][i]);
1184 }
1185 return trac;
1186 }
1187
1188 // DETERMINANT
1189 // Returns the determinant of a complex square matrix [instance method]
1190 public Complex determinant(){
1191 int n = this.nrow;
1192 if(n!=this.ncol)throw new IllegalArgumentException("Matrix is not square");
1193 Complex det = new Complex();
1194 ComplexMatrix ludmat;
1195
1196 ludmat = this.luDecomp();
1197 det.reset(ludmat.dswap,0.0);
1198 for(int j=0; j<n; j++){
1199 det.timesEquals(ludmat.matrix[j][j]);
1200 }
1201 return det;
1202 }
1203
1204 // Returns the determinant of a complex square matrix [static method]
1205 public static Complex determinant(ComplexMatrix amat){
1206 int n = amat.nrow;
1207 if(n!=amat.ncol)throw new IllegalArgumentException("Matrix is not square");
1208 Complex det = new Complex();
1209 ComplexMatrix ludmat;
1210
1211 ludmat = amat.luDecomp();
1212 det.reset(ludmat.dswap,0.0);
1213 for(int j=0; j<n; j++){
1214 det.timesEquals(ludmat.matrix[j][j]);
1215 }
1216 return det;
1217 }
1218
1219 // Returns the log(determinant) of a complex square matrix [instance method].
1220 // Useful if determinant() underflows or overflows.
1221 public Complex logDeterminant(){
1222 int n = this.nrow;
1223 if(n!=this.ncol)throw new IllegalArgumentException("Matrix is not square");
1224 Complex det = new Complex();
1225 ComplexMatrix ludmat;
1226
1227 ludmat = this.luDecomp();
1228 det.reset(ludmat.dswap,0.0);
1229 det=Complex.log(det);
1230 for(int j=0; j<n; j++){
1231 det.plusEquals(Complex.log(ludmat.matrix[j][j]));
1232 }
1233 return det;
1234 }
1235
1236 // Returns the log(determinant) of a complex square matrix [static method].
1237 // Useful if determinant() underflows or overflows.
1238 public static Complex logDeterminant(ComplexMatrix amat){
1239 int n = amat.nrow;
1240 if(n!=amat.ncol)throw new IllegalArgumentException("Matrix is not square");
1241 Complex det = new Complex();
1242 ComplexMatrix ludmat;
1243
1244 ludmat = amat.luDecomp();
1245 det.reset(ludmat.dswap,0.0);
1246 det=Complex.log(det);
1247 for(int j=0; j<n; j++){
1248 det.plusEquals(Complex.log(ludmat.matrix[j][j]));
1249 }
1250 return det;
1251 }
1252
1253 // REDUCED ROW ECHELON FORM
1254 public ComplexMatrix reducedRowEchelonForm() {
1255
1256 Complex[][] mat = Complex.twoDarray(this.nrow, this.ncol);
1257 for(int i=0; i<this.nrow; i++){
1258 for(int j=0; j<this.ncol; j++){
1259 mat[i][j] = this.matrix[i][j];
1260 }
1261 }
1262
1263 int leadingCoeff = 0;
1264 int rowPointer = 0;
1265
1266 boolean testOuter = true;
1267 while(testOuter){
1268 int counter = rowPointer;
1269 boolean testInner = true;
1270 while(testInner && mat[counter][leadingCoeff].equals(Complex.zero())) {
1271 counter++;
1272 if(counter == this.nrow){
1273 counter = rowPointer;
1274 leadingCoeff++;
1275 if(leadingCoeff == this.ncol)testInner=false;
1276 }
1277 }
1278 if(testInner){
1279 Complex[] temp = mat[rowPointer];
1280 mat[rowPointer] = mat[counter];
1281 mat[counter] = temp;
1282
1283 Complex pivot = mat[rowPointer][leadingCoeff];
1284 for(int j=0; j<this.ncol; j++)mat[rowPointer][j] = mat[rowPointer][j].over(pivot);
1285
1286 for(int i=0; i<this.nrow; i++){
1287 if (i!=rowPointer) {
1288 pivot = mat[i][leadingCoeff];
1289 for (int j=0; j<this.ncol; j++)mat[i][j] = mat[i][j].minus(pivot.times(mat[rowPointer][j]));
1290 }
1291 }
1292 leadingCoeff++;
1293 if(leadingCoeff>=this.ncol)testOuter = false;
1294 }
1295 rowPointer++;
1296 if(rowPointer>=this.nrow || !testInner)testOuter = false;
1297 }
1298
1299 for(int i=0; i<this.nrow; i++){
1300 for(int j=0; j<this.ncol; j++){
1301 if(mat[i][j].getReal()==-0.0)mat[i][j].reset(0.0, mat[i][j].getImag());
1302 if(mat[i][j].getImag()==-0.0)mat[i][j].reset(mat[i][j].getReal(), 0.0);
1303 }
1304 }
1305 return new ComplexMatrix(mat);
1306 }
1307
1308 // FROBENIUS NORM of a complex matrix
1309 // Sometimes referred to as the EUCLIDEAN NORM
1310 public double frobeniusNorm(){
1311 double norm=0.0D;
1312 for(int i=0; i<this.nrow; i++){
1313 for(int j=0; j<this.ncol; j++){
1314 norm=Fmath.hypot(norm, Complex.abs(matrix[i][j]));
1315 }
1316 }
1317 return norm;
1318 }
1319
1320 // ONE NORM of a complex matrix
1321 public double oneNorm(){
1322 double norm=0.0D;
1323 double sum = 0.0D;
1324 for(int i=0; i<this.nrow; i++){
1325 sum=0.0D;
1326 for(int j=0; j<this.ncol; j++){
1327 sum+=Complex.abs(this.matrix[i][j]);
1328 }
1329 norm=Math.max(norm,sum);
1330 }
1331 return norm;
1332 }
1333
1334 // INFINITY NORM of a complex matrix
1335 public double infinityNorm(){
1336 double norm=0.0D;
1337 double sum=0.0D;
1338 for(int i=0; i<this.nrow; i++){
1339 sum=0.0D;
1340 for(int j=0; j<this.ncol; j++){
1341 sum+=Complex.abs(this.matrix[i][j]);
1342 }
1343 norm=Math.max(norm,sum);
1344 }
1345 return norm;
1346 }
1347
1348
1349 // LU DECOMPOSITION OF COMPLEX MATRIX A
1350 // For details of LU decomposition
1351 // See Numerical Recipes, The Art of Scientific Computing
1352 // by W H Press, S A Teukolsky, W T Vetterling & B P Flannery
1353 // Cambridge University Press, http://www.nr.com/
1354 // ComplexMatrix ludmat is the returned LU decompostion
1355 // int[] index is the vector of row permutations
1356 // dswap returns +1.0 for even number of row interchanges
1357 // returns -1.0 for odd number of row interchanges
1358 public ComplexMatrix luDecomp(){
1359 if(this.nrow!=this.ncol)throw new IllegalArgumentException("A matrix is not square");
1360 int n=this.nrow;
1361 int imax=0;
1362 double dum=0.0D, temp=0.0D, big=0.0D;
1363 double[] vv = new double[n];
1364 Complex sum = new Complex();
1365 Complex dumm = new Complex();
1366
1367 ComplexMatrix ludmat=ComplexMatrix.copy(this);
1368 Complex[][] ludarray = ludmat.getArrayReference();
1369
1370 ludmat.dswap=1.0;
1371 for (int i=0;i<n;i++) {
1372 big=0.0;
1373 for (int j=0;j<n;j++){
1374 if ((temp=Complex.abs(ludarray[i][j])) > big) big=temp;
1375 }
1376 if (big == 0.0) throw new ArithmeticException("Singular matrix");
1377 vv[i]=1.0/big;
1378 }
1379 for (int j=0;j<n;j++) {
1380 for (int i=0;i<j;i++) {
1381 sum=Complex.copy(ludarray[i][j]);
1382 for (int k=0;k<i;k++) sum.minusEquals(ludarray[i][k].times(ludarray[k][j]));
1383 ludarray[i][j]=Complex.copy(sum);
1384 }
1385 big=0.0;
1386 for (int i=j;i<n;i++) {
1387 sum=Complex.copy(ludarray[i][j]);
1388 for (int k=0;k<j;k++){
1389 sum.minusEquals(ludarray[i][k].times(ludarray[k][j]));
1390 }
1391 ludarray[i][j]=Complex.copy(sum);
1392 if ((dum=vv[i]*Complex.abs(sum)) >= big) {
1393 big=dum;
1394 imax=i;
1395 }
1396 }
1397 if (j != imax) {
1398 for (int k=0;k<n;k++) {
1399 dumm=Complex.copy(ludarray[imax][k]);
1400 ludarray[imax][k]=Complex.copy(ludarray[j][k]);
1401 ludarray[j][k]=Complex.copy(dumm);
1402 }
1403 ludmat.dswap = -ludmat.dswap;
1404 vv[imax]=vv[j];
1405 }
1406 ludmat.index[j]=imax;
1407
1408 if(ludarray[j][j].isZero()){
1409 ludarray[j][j].reset(TINY, TINY);
1410 }
1411 if(j != n-1) {
1412 dumm=Complex.over(1.0,ludarray[j][j]);
1413 for (int i=j+1;i<n;i++){
1414 ludarray[i][j].timesEquals(dumm);
1415 }
1416 }
1417 }
1418 return ludmat;
1419 }
1420
1421 // Solves the set of n linear complex equations A.X=B using not A but its LU decomposition
1422 // Complex bvec is the vector B (input)
1423 // Complex xvec is the vector X (output)
1424 // index is the permutation vector produced by luDecomp()
1425 public Complex[] luBackSub(Complex[] bvec){
1426 int ii=0,ip=0;
1427 int n=bvec.length;
1428 if(n!=this.ncol)throw new IllegalArgumentException("vector length is not equal to matrix dimension");
1429 if(this.ncol!=this.nrow)throw new IllegalArgumentException("matrix is not square");
1430 Complex sum=new Complex();
1431 Complex[] xvec=new Complex[n];
1432 for(int i=0; i<n; i++){
1433 xvec[i]=Complex.copy(bvec[i]);
1434 }
1435 for (int i=0;i<n;i++) {
1436 ip=this.index[i];
1437 sum=Complex.copy(xvec[ip]);
1438 xvec[ip]=Complex.copy(xvec[i]);
1439 if (ii==0){
1440 for (int j=ii;j<=i-1;j++){
1441 sum.minusEquals(this.matrix[i][j].times(xvec[j]));
1442 }
1443 }
1444 else{
1445 if(sum.isZero()) ii=i;
1446 }
1447 xvec[i]=Complex.copy(sum);
1448 }
1449 for(int i=n-1;i>=0;i--) {
1450 sum=Complex.copy(xvec[i]);
1451 for (int j=i+1;j<n;j++){
1452 sum.minusEquals(this.matrix[i][j].times(xvec[j]));
1453 }
1454 xvec[i]= sum.over(this.matrix[i][i]);
1455 }
1456 return xvec;
1457 }
1458
1459 // Solves the set of n linear complex equations A.X=B
1460 // Complex bvec is the vector B (input)
1461 // Complex xvec is the vector X (output)
1462 public Complex[] solveLinearSet(Complex[] bvec){
1463 ComplexMatrix ludmat;
1464
1465 ludmat=this.luDecomp();
1466 return ludmat.luBackSub(bvec);
1467 }
1468}
1469
Note: See TracBrowser for help on using the repository browser.