|
1
2 package JSci.maths.matrices;
3
4 import JSci.maths.ArrayMath;
5 import JSci.maths.Complex;
6 import JSci.maths.ComplexMapping;
7 import JSci.maths.LinearMath;
8 import JSci.maths.DimensionException;
9 import JSci.maths.MaximumIterationsExceededException;
10 import JSci.maths.vectors.AbstractComplexVector;
11 import JSci.maths.vectors.ComplexVector;
12
13
19 public class ComplexTridiagonalMatrix extends AbstractComplexSquareMatrix implements TridiagonalMatrix {
20
23 protected final double ldiagRe[],ldiagIm[];
24 protected final double diagRe[],diagIm[];
25 protected final double udiagRe[],udiagIm[];
26
30 public ComplexTridiagonalMatrix(final int size) {
31 super(size);
32 ldiagRe=new double[size];
33 ldiagIm=new double[size];
34 diagRe=new double[size];
35 diagIm=new double[size];
36 udiagRe=new double[size];
37 udiagIm=new double[size];
38 }
39
45 public ComplexTridiagonalMatrix(final Complex array[][]) {
46 this(array.length);
47 if(!ArrayMath.isSquare(array))
48 throw new MatrixDimensionException("Array is not square.");
49 diagRe[0]=array[0][0].real();
50 diagIm[0]=array[0][0].imag();
51 udiagRe[0]=array[0][1].real();
52 udiagIm[0]=array[0][1].imag();
53 int i=1;
54 for(;i<array.length-1;i++) {
55 ldiagRe[i]=array[i][i-1].real();
56 ldiagIm[i]=array[i][i-1].imag();
57 diagRe[i]=array[i][i].real();
58 diagIm[i]=array[i][i].imag();
59 udiagRe[i]=array[i][i+1].real();
60 udiagIm[i]=array[i][i+1].imag();
61 }
62 ldiagRe[i]=array[i][i-1].real();
63 ldiagIm[i]=array[i][i-1].imag();
64 diagRe[i]=array[i][i].real();
65 diagIm[i]=array[i][i].imag();
66 }
67
71 public boolean equals(AbstractComplexMatrix m, double tol) {
72 if(m instanceof TridiagonalMatrix) {
73 if(numRows != m.rows() || numCols != m.columns())
74 return false;
75 double sumSqr = 0;
76 double ldeltaRe,ldeltaIm;
77 double deltaRe = diagRe[0] - m.getRealElement(0,0);
78 double deltaIm = diagIm[0] - m.getImagElement(0,0);
79 double udeltaRe = udiagRe[0] - m.getRealElement(0,1);
80 double udeltaIm = udiagIm[0] - m.getImagElement(0,1);
81 sumSqr += deltaRe*deltaRe+deltaIm*deltaIm + udeltaRe*udeltaRe+udeltaIm*udeltaIm;
82 int i=1;
83 for(;i<numRows-1;i++) {
84 ldeltaRe = ldiagRe[i] - m.getRealElement(i,i-1);
85 ldeltaIm = ldiagIm[i] - m.getImagElement(i,i-1);
86 deltaRe = diagRe[i] - m.getRealElement(i,i);
87 deltaIm = diagIm[i] - m.getImagElement(i,i);
88 udeltaRe = udiagRe[i] - m.getRealElement(i,i+1);
89 udeltaIm = udiagIm[i] - m.getImagElement(i,i+1);
90 sumSqr += ldeltaRe*ldeltaRe+ldeltaIm*ldeltaIm + deltaRe*deltaRe+deltaIm*deltaIm + udeltaRe*udeltaRe+udeltaIm*udeltaIm;
91 }
92 ldeltaRe = ldiagRe[i] - m.getRealElement(i,i-1);
93 ldeltaIm = ldiagIm[i] - m.getImagElement(i,i-1);
94 deltaRe = diagRe[i] - m.getRealElement(i,i);
95 deltaIm = diagIm[i] - m.getImagElement(i,i);
96 sumSqr += ldeltaRe*ldeltaRe+ldeltaIm*ldeltaIm + deltaRe*deltaRe+deltaIm*deltaIm;
97 return (sumSqr <= tol*tol);
98 } else {
99 return false;
100 }
101 }
102
105 public String  toString() {
106 final StringBuffer buf=new StringBuffer (5*rows()*columns());
107 for(int i=0;i<rows();i++) {
108 for(int j=0;j<columns();j++) {
109 buf.append(getElement(i,j).toString());
110 buf.append(' ');
111 }
112 buf.append('\n');
113 }
114 return buf.toString();
115 }
116
120 public AbstractDoubleMatrix real() {
121 final DoubleTridiagonalMatrix m=new DoubleTridiagonalMatrix(numRows);
122 m.diag[0]=diagRe[0];
123 m.udiag[0]=udiagRe[0];
124 int i=1;
125 for(;i<numRows-1;i++) {
126 m.ldiag[i]=ldiagRe[i];
127 m.diag[i]=diagRe[i];
128 m.udiag[i]=udiagRe[i];
129 }
130 m.ldiag[i]=ldiagRe[i];
131 m.diag[i]=diagRe[i];
132 return m;
133 }
134
138 public AbstractDoubleMatrix imag() {
139 final DoubleTridiagonalMatrix m=new DoubleTridiagonalMatrix(numRows);
140 m.diag[0]=diagIm[0];
141 m.udiag[0]=udiagIm[0];
142 int i=1;
143 for(;i<numRows-1;i++) {
144 m.ldiag[i]=ldiagIm[i];
145 m.diag[i]=diagIm[i];
146 m.udiag[i]=udiagIm[i];
147 }
148 m.ldiag[i]=ldiagIm[i];
149 m.diag[i]=diagIm[i];
150 return m;
151 }
152
158 public Complex getElement(final int i, final int j) {
159 if(i>=0 && i<numRows && j>=0 && j<numCols) {
160 if(j == i-1)
161 return new Complex(ldiagRe[i], ldiagIm[i]);
162 else if(j == i)
163 return new Complex(diagRe[i], diagIm[i]);
164 else if(j == i+1)
165 return new Complex(udiagRe[i], udiagIm[i]);
166 else
167 return Complex.ZERO;
168 } else
169 throw new MatrixDimensionException(getInvalidElementMsg(i,j));
170 }
171 public double getRealElement(final int i, final int j) {
172 if(i>=0 && i<numRows && j>=0 && j<numCols) {
173 if(j == i-1)
174 return ldiagRe[i];
175 else if(j == i)
176 return diagRe[i];
177 else if(j == i+1)
178 return udiagRe[i];
179 else
180 return 0.0;
181 } else
182 throw new MatrixDimensionException(getInvalidElementMsg(i,j));
183 }
184 public double getImagElement(final int i, final int j) {
185 if(i>=0 && i<numRows && j>=0 && j<numCols) {
186 if(j == i-1)
187 return ldiagIm[i];
188 else if(j == i)
189 return diagIm[i];
190 else if(j == i+1)
191 return udiagIm[i];
192 else
193 return 0.0;
194 } else
195 throw new MatrixDimensionException(getInvalidElementMsg(i,j));
196 }
197
205 public void setElement(final int i, final int j, final Complex z) {
206 if(i>=0 && i<numRows && j>=0 && j<numCols) {
207 if(j == i-1) {
208 ldiagRe[i] = z.real();
209 ldiagIm[i] = z.imag();
210 } else if(j == i) {
211 diagRe[i] = z.real();
212 diagIm[i] = z.imag();
213 } else if(j == i+1) {
214 udiagRe[i] = z.real();
215 udiagIm[i] = z.imag();
216 } else {
217 throw new MatrixDimensionException(getInvalidElementMsg(i,j));
218 }
219 } else
220 throw new MatrixDimensionException(getInvalidElementMsg(i,j));
221 }
222
231 public void setElement(final int i, final int j, final double x, final double y) {
232 if(i>=0 && i<numRows && j>=0 && j<numCols) {
233 if(j == i-1) {
234 ldiagRe[i] = x;
235 ldiagIm[i] = y;
236 } else if(j == i) {
237 diagRe[i] = x;
238 diagIm[i] = y;
239 } else if(j == i+1) {
240 udiagRe[i] = x;
241 udiagIm[i] = y;
242 } else {
243 throw new MatrixDimensionException(getInvalidElementMsg(i,j));
244 }
245 } else
246 throw new MatrixDimensionException(getInvalidElementMsg(i,j));
247 }
248
251 public Complex trace() {
252 double trRe=diagRe[0];
253 double trIm=diagIm[0];
254 for(int i=1;i<numRows;i++) {
255 trRe+=diagRe[i];
256 trIm+=diagIm[i];
257 }
258 return new Complex(trRe,trIm);
259 }
260
264 public double infNorm() {
265 double result = Math.sqrt((diagRe[0]*diagRe[0] + diagIm[0]*diagIm[0])) + Math.sqrt((udiagRe[0]*udiagRe[0] + udiagIm[0]*udiagIm[0]));
266 double tmpResult;
267 int i=1;
268 for(;i<numRows-1;i++) {
269 tmpResult = Math.sqrt((ldiagRe[i]*ldiagRe[i] + ldiagIm[i]*ldiagIm[i])) + Math.sqrt((diagRe[i]*diagRe[i] + diagIm[i]*diagIm[i])) + Math.sqrt((udiagRe[i]*udiagRe[i] + udiagIm[i]*udiagIm[i]));
270 if(tmpResult>result)
271 result=tmpResult;
272 }
273 tmpResult = Math.sqrt((ldiagRe[i]*ldiagRe[i] + ldiagIm[i]*ldiagIm[i])) + Math.sqrt((diagRe[i]*diagRe[i] + diagIm[i]*diagIm[i]));
274 if(tmpResult>result)
275 result=tmpResult;
276 return result;
277 }
278
282 public double frobeniusNorm() {
283 double result=diagRe[0]*diagRe[0]+diagIm[0]*diagIm[0]+
284 udiagRe[0]*udiagRe[0]+udiagIm[0]*udiagIm[0];
285 int i=1;
286 for(;i<numRows-1;i++) {
287 result+=ldiagRe[i]*ldiagRe[i]+ldiagIm[i]*ldiagIm[i]+
288 diagRe[i]*diagRe[i]+diagIm[i]*diagIm[i]+
289 udiagRe[i]*udiagRe[i]+udiagIm[i]*udiagIm[i];
290 }
291 result+=ldiagRe[i]*ldiagRe[i]+ldiagIm[i]*ldiagIm[i]+
292 diagRe[i]*diagRe[i]+diagIm[i]*diagIm[i];
293 return Math.sqrt(result);
294 }
295
299 public double operatorNorm() throws MaximumIterationsExceededException {
300 return Math.sqrt(ArrayMath.max(LinearMath.eigenvalueSolveHermitian((ComplexTridiagonalMatrix)(this.hermitianAdjoint().multiply(this)))));
301 }
302
303
307
309
314 public AbstractComplexSquareMatrix add(final AbstractComplexSquareMatrix m) {
315 if(m instanceof ComplexTridiagonalMatrix)
316 return add((ComplexTridiagonalMatrix)m);
317 if(m instanceof TridiagonalMatrix)
318 return addTridiagonal(m);
319 if(m instanceof ComplexSquareMatrix)
320 return add((ComplexSquareMatrix)m);
321
322 if(numRows==m.rows() && numCols==m.columns()) {
323 final double arrayRe[][]=new double[numRows][numCols];
324 final double arrayIm[][]=new double[numRows][numCols];
325 for(int i=0;i<numRows;i++) {
326 Complex elem=getElement(i,0).add(m.getElement(i,0));
327 arrayRe[i][0]=elem.real();
328 arrayIm[i][0]=elem.imag();
329 for(int j=1;j<numCols;j++) {
330 elem=getElement(i,j).add(m.getElement(i,j));
331 arrayRe[i][j]=elem.real();
332 arrayIm[i][j]=elem.imag();
333 }
334 }
335 return new ComplexSquareMatrix(arrayRe,arrayIm);
336 } else {
337 throw new MatrixDimensionException("Matrices are different sizes.");
338 }
339 }
340 public ComplexSquareMatrix add(final ComplexSquareMatrix m) {
341 if(numRows==m.numRows && numCols==m.numCols) {
342 final double arrayRe[][]=new double[numRows][numRows];
343 final double arrayIm[][]=new double[numRows][numRows];
344 for(int i=0;i<numRows;i++) {
345 System.arraycopy(m.matrixRe[i],0,arrayRe[i],0,numRows);
346 System.arraycopy(m.matrixIm[i],0,arrayIm[i],0,numRows);
347 }
348 arrayRe[0][0]+=diagRe[0];
349 arrayIm[0][0]+=diagIm[0];
350 arrayRe[0][1]+=udiagRe[0];
351 arrayIm[0][1]+=udiagIm[0];
352 int n=numCols-1;
353 for(int i=1;i<n;i++) {
354 arrayRe[i][i-1]+=ldiagRe[i];
355 arrayIm[i][i-1]+=ldiagIm[i];
356 arrayRe[i][i]+=diagRe[i];
357 arrayIm[i][i]+=diagIm[i];
358 arrayRe[i][i+1]+=udiagRe[i];
359 arrayIm[i][i+1]+=udiagIm[i];
360 }
361 arrayRe[n][n-1]+=ldiagRe[n];
362 arrayIm[n][n-1]+=ldiagIm[n];
363 arrayRe[n][n]+=diagRe[n];
364 arrayIm[n][n]+=diagIm[n];
365 return new ComplexSquareMatrix(arrayRe,arrayIm);
366 } else
367 throw new MatrixDimensionException("Matrices are different sizes.");
368 }
369 public ComplexTridiagonalMatrix add(final ComplexTridiagonalMatrix m) {
370 int mRow=numRows;
371 if(mRow==m.numRows) {
372 final ComplexTridiagonalMatrix ans=new ComplexTridiagonalMatrix(mRow);
373 ans.diagRe[0]=diagRe[0]+m.diagRe[0];
374 ans.diagIm[0]=diagIm[0]+m.diagIm[0];
375 ans.udiagRe[0]=udiagRe[0]+m.udiagRe[0];
376 ans.udiagIm[0]=udiagIm[0]+m.udiagIm[0];
377 mRow--;
378 for(int i=1;i<mRow;i++) {
379 ans.ldiagRe[i]=ldiagRe[i]+m.ldiagRe[i];
380 ans.ldiagIm[i]=ldiagIm[i]+m.ldiagIm[i];
381 ans.diagRe[i]=diagRe[i]+m.diagRe[i];
382 ans.diagIm[i]=diagIm[i]+m.diagIm[i];
383 ans.udiagRe[i]=udiagRe[i]+m.udiagRe[i];
384 ans.udiagIm[i]=udiagIm[i]+m.udiagIm[i];
385 }
386 ans.ldiagRe[mRow]=ldiagRe[mRow]+m.ldiagRe[mRow];
387 ans.ldiagIm[mRow]=ldiagIm[mRow]+m.ldiagIm[mRow];
388 ans.diagRe[mRow]=diagRe[mRow]+m.diagRe[mRow];
389 ans.diagIm[mRow]=diagIm[mRow]+m.diagIm[mRow];
390 return ans;
391 } else
392 throw new MatrixDimensionException("Matrices are different sizes.");
393 }
394 private ComplexTridiagonalMatrix addTridiagonal(final AbstractComplexSquareMatrix m) {
395 int mRow=numRows;
396 if(mRow==m.rows()) {
397 final ComplexTridiagonalMatrix ans=new ComplexTridiagonalMatrix(mRow);
398 Complex elem=m.getElement(0,0);
399 ans.diagRe[0]=diagRe[0]+elem.real();
400 ans.diagIm[0]=diagIm[0]+elem.imag();
401 elem=m.getElement(0,1);
402 ans.udiagRe[0]=udiagRe[0]+elem.real();
403 ans.udiagIm[0]=udiagIm[0]+elem.imag();
404 mRow--;
405 for(int i=1;i<mRow;i++) {
406 elem=m.getElement(i,i-1);
407 ans.ldiagRe[i]=ldiagRe[i]+elem.real();
408 ans.ldiagIm[i]=ldiagIm[i]+elem.imag();
409 elem=m.getElement(i,i);
410 ans.diagRe[i]=diagRe[i]+elem.real();
411 ans.diagIm[i]=diagIm[i]+elem.imag();
412 elem=m.getElement(i,i+1);
413 ans.udiagRe[i]=udiagRe[i]+elem.real();
414 ans.udiagIm[i]=udiagIm[i]+elem.imag();
415 }
416 elem=m.getElement(mRow,mRow-1);
417 ans.ldiagRe[mRow]=ldiagRe[mRow]+elem.real();
418 ans.ldiagIm[mRow]=ldiagIm[mRow]+elem.imag();
419 elem=m.getElement(mRow,mRow);
420 ans.diagRe[mRow]=diagRe[mRow]+elem.real();
421 ans.diagIm[mRow]=diagIm[mRow]+elem.imag();
422 return ans;
423 } else {
424 throw new MatrixDimensionException("Matrices are different sizes.");
425 }
426 }
427
428
430
435 public AbstractComplexSquareMatrix subtract(final AbstractComplexSquareMatrix m) {
436 if(m instanceof ComplexTridiagonalMatrix)
437 return subtract((ComplexTridiagonalMatrix)m);
438 if(m instanceof TridiagonalMatrix)
439 return subtractTridiagonal(m);
440 if(m instanceof ComplexSquareMatrix)
441 return subtract((ComplexSquareMatrix)m);
442
443 if(numRows==m.rows() && numCols==m.columns()) {
444 final double arrayRe[][]=new double[numRows][numCols];
445 final double arrayIm[][]=new double[numRows][numCols];
446 Complex elem;
447 for(int i=0;i<numRows;i++) {
448 elem=getElement(i,0).subtract(m.getElement(i,0));
449 arrayRe[i][0]=elem.real();
450 arrayIm[i][0]=elem.imag();
451 for(int j=1;j<numCols;j++) {
452 elem=getElement(i,j).subtract(m.getElement(i,j));
453 arrayRe[i][j]=elem.real();
454 arrayIm[i][j]=elem.imag();
455 }
456 }
457 return new ComplexSquareMatrix(arrayRe,arrayIm);
458 } else {
459 throw new MatrixDimensionException("Matrices are different sizes.");
460 }
461 }
462 public ComplexSquareMatrix subtract(final ComplexSquareMatrix m) {
463 if(numRows==m.numRows && numCols==m.numCols) {
464 final double arrayRe[][]=new double[numRows][numCols];
465 final double arrayIm[][]=new double[numRows][numCols];
466 Complex elem;
467 for(int j,i=0;i<numRows;i++) {
468 elem=getElement(i,0);
469 arrayRe[i][0]=elem.real()-m.matrixRe[i][0];
470 arrayIm[i][0]=elem.imag()-m.matrixIm[i][0];
471 for(j=1;j<numCols;j++) {
472 elem=getElement(i,j);
473 arrayRe[i][j]=elem.real()-m.matrixRe[i][j];
474 arrayIm[i][j]=elem.imag()-m.matrixIm[i][j];
475 }
476 }
477 return new ComplexSquareMatrix(arrayRe,arrayIm);
478 } else
479 throw new MatrixDimensionException("Matrices are different sizes.");
480 }
481
486 public ComplexTridiagonalMatrix subtract(final ComplexTridiagonalMatrix m) {
487 int mRow=numRows;
488 if(mRow==m.numRows) {
489 final ComplexTridiagonalMatrix ans=new ComplexTridiagonalMatrix(mRow);
490 ans.diagRe[0]=diagRe[0]-m.diagRe[0];
491 ans.diagIm[0]=diagIm[0]-m.diagIm[0];
492 ans.udiagRe[0]=udiagRe[0]-m.udiagRe[0];
493 ans.udiagIm[0]=udiagIm[0]-m.udiagIm[0];
494 mRow--;
495 for(int i=1;i<mRow;i++) {
496 ans.ldiagRe[i]=ldiagRe[i]-m.ldiagRe[i];
497 ans.ldiagIm[i]=ldiagIm[i]-m.ldiagIm[i];
498 ans.diagRe[i]=diagRe[i]-m.diagRe[i];
499 ans.diagIm[i]=diagIm[i]-m.diagIm[i];
500 ans.udiagRe[i]=udiagRe[i]-m.udiagRe[i];
501 ans.udiagIm[i]=udiagIm[i]-m.udiagIm[i];
502 }
503 ans.ldiagRe[mRow]=ldiagRe[mRow]-m.ldiagRe[mRow];
504 ans.ldiagIm[mRow]=ldiagIm[mRow]-m.ldiagIm[mRow];
505 ans.diagRe[mRow]=diagRe[mRow]-m.diagRe[mRow];
506 ans.diagIm[mRow]=diagIm[mRow]-m.diagIm[mRow];
507 return ans;
508 } else
509 throw new MatrixDimensionException("Matrices are different sizes.");
510 }
511 private ComplexTridiagonalMatrix subtractTridiagonal(final AbstractComplexSquareMatrix m) {
512 int mRow=numRows;
513 if(mRow==m.rows()) {
514 final ComplexTridiagonalMatrix ans=new ComplexTridiagonalMatrix(mRow);
515 Complex elem=m.getElement(0,0);
516 ans.diagRe[0]=diagRe[0]-elem.real();
517 ans.diagIm[0]=diagIm[0]-elem.imag();
518 elem=m.getElement(0,1);
519 ans.udiagRe[0]=udiagRe[0]-elem.real();
520 ans.udiagIm[0]=udiagIm[0]-elem.imag();
521 mRow--;
522 for(int i=1;i<mRow;i++) {
523 elem=m.getElement(i,i-1);
524 ans.ldiagRe[i]=ldiagRe[i]-elem.real();
525 ans.ldiagIm[i]=ldiagIm[i]-elem.imag();
526 elem=m.getElement(i,i);
527 ans.diagRe[i]=diagRe[i]-elem.real();
528 ans.diagIm[i]=diagIm[i]-elem.imag();
529 elem=m.getElement(i,i+1);
530 ans.udiagRe[i]=udiagRe[i]-elem.real();
531 ans.udiagIm[i]=udiagIm[i]-elem.imag();
532 }
533 elem=m.getElement(mRow,mRow-1);
534 ans.ldiagRe[mRow]=ldiagRe[mRow]-elem.real();
535 ans.ldiagIm[mRow]=ldiagIm[mRow]-elem.imag();
536 elem=m.getElement(mRow,mRow);
537 ans.diagRe[mRow]=diagRe[mRow]-elem.real();
538 ans.diagIm[mRow]=diagIm[mRow]-elem.imag();
539 return ans;
540 } else {
541 throw new MatrixDimensionException("Matrices are different sizes.");
542 }
543 }
544
545
547
552 public AbstractComplexMatrix scalarMultiply(final Complex z) {
553 final double real=z.real();
554 final double imag=z.imag();
555 int mRow=numRows;
556 final ComplexTridiagonalMatrix ans=new ComplexTridiagonalMatrix(mRow);
557 ans.diagRe[0]=real*diagRe[0]-imag*diagIm[0];
558 ans.diagIm[0]=imag*diagRe[0]+real*diagIm[0];
559 ans.udiagRe[0]=real*udiagRe[0]-imag*udiagIm[0];
560 ans.udiagIm[0]=imag*udiagRe[0]+real*udiagIm[0];
561 mRow--;
562 for(int i=1;i<mRow;i++) {
563 ans.ldiagRe[i]=real*ldiagRe[i]-imag*ldiagIm[i];
564 ans.ldiagIm[i]=imag*ldiagRe[i]+real*ldiagIm[i];
565 ans.diagRe[i]=real*diagRe[i]-imag*diagIm[i];
566 ans.diagIm[i]=imag*diagRe[i]+real*diagIm[i];
567 ans.udiagRe[i]=real*udiagRe[i]-imag*udiagIm[i];
568 ans.udiagIm[i]=imag*udiagRe[i]+real*udiagIm[i];
569 }
570 ans.ldiagRe[mRow]=real*ldiagRe[mRow]-imag*ldiagIm[mRow];
571 ans.ldiagIm[mRow]=imag*ldiagRe[mRow]+real*ldiagIm[mRow];
572 ans.diagRe[mRow]=real*diagRe[mRow]-imag*diagIm[mRow];
573 ans.diagIm[mRow]=imag*diagRe[mRow]+real*diagIm[mRow];
574 return ans;
575 }
576
581 public AbstractComplexMatrix scalarMultiply(final double x) {
582 int mRow=numRows;
583 final ComplexTridiagonalMatrix ans=new ComplexTridiagonalMatrix(mRow);
584 ans.diagRe[0]=x*diagRe[0];
585 ans.diagIm[0]=x*diagIm[0];
586 ans.udiagRe[0]=x*udiagRe[0];
587 ans.udiagIm[0]=x*udiagIm[0];
588 mRow--;
589 for(int i=1;i<mRow;i++) {
590 ans.ldiagRe[i]=x*ldiagRe[i];
591 ans.ldiagIm[i]=x*ldiagIm[i];
592 ans.diagRe[i]=x*diagRe[i];
593 ans.diagIm[i]=x*diagIm[i];
594 ans.udiagRe[i]=x*udiagRe[i];
595 ans.udiagIm[i]=x*udiagIm[i];
596 }
597 ans.ldiagRe[mRow]=x*ldiagRe[mRow];
598 ans.ldiagIm[mRow]=x*ldiagIm[mRow];
599 ans.diagRe[mRow]=x*diagRe[mRow];
600 ans.diagIm[mRow]=x*diagIm[mRow];
601 return ans;
602 }
603
604
606
611 public AbstractComplexVector multiply(final AbstractComplexVector v) {
612 int mRow=numRows;
613 if(mRow==v.dimension()) {
614 final double arrayRe[]=new double[mRow];
615 final double arrayIm[]=new double[mRow];
616 Complex comp;
617 comp=v.getComponent(0);
618 arrayRe[0]=(diagRe[0]*comp.real() - diagIm[0]*comp.imag());
619 arrayIm[0]=(diagIm[0]*comp.real() + diagRe[0]*comp.imag());
620 comp=v.getComponent(1);
621 arrayRe[0]+=(udiagRe[0]*comp.real() - udiagIm[0]*comp.imag());
622 arrayIm[0]+=(udiagIm[0]*comp.real() + udiagRe[0]*comp.imag());
623 mRow--;
624 for(int i=1;i<mRow;i++) {
625 comp=v.getComponent(i-1);
626 arrayRe[i]=(ldiagRe[i]*comp.real() - ldiagIm[i]*comp.imag());
627 arrayIm[i]=(ldiagIm[i]*comp.real() + ldiagRe[i]*comp.imag());
628 comp=v.getComponent(i);
629 arrayRe[i]+=(diagRe[i]*comp.real() - diagIm[i]*comp.imag());
630 arrayIm[i]+=(diagIm[i]*comp.real() + diagRe[i]*comp.imag());
631 comp=v.getComponent(i+1);
632 arrayRe[i]+=(udiagRe[i]*comp.real() - udiagIm[i]*comp.imag());
633 arrayIm[i]+=(udiagIm[i]*comp.real() + udiagRe[i]*comp.imag());
634 }
635 comp=v.getComponent(mRow-1);
636 arrayRe[mRow]=(ldiagRe[mRow]*comp.real() - ldiagIm[mRow]*comp.imag());
637 arrayIm[mRow]=(ldiagIm[mRow]*comp.real() + ldiagRe[mRow]*comp.imag());
638 comp=v.getComponent(mRow);
639 arrayRe[mRow]+=(diagRe[mRow]*comp.real() - diagIm[mRow]*comp.imag());
640 arrayIm[mRow]+=(diagIm[mRow]*comp.real() + diagRe[mRow]*comp.imag());
641 return new ComplexVector(arrayRe,arrayIm);
642 } else
643 throw new DimensionException("Matrix and vector are incompatible.");
644 }
645
650 public AbstractComplexSquareMatrix multiply(final AbstractComplexSquareMatrix m) {
651 if(m instanceof ComplexTridiagonalMatrix)
652 return multiply((ComplexTridiagonalMatrix)m);
653 if(m instanceof TridiagonalMatrix)
654 return multiplyTridiagonal(m);
655 if(m instanceof ComplexSquareMatrix)
656 return multiply((ComplexSquareMatrix)m);
657
658 if(numCols==m.rows()) {
659 final int mColumns = m.columns();
660 final double arrayRe[][]=new double[numRows][mColumns];
661 final double arrayIm[][]=new double[numRows][mColumns];
662 final int lastRow = numRows-1;
663 for(int j=0; j<numRows; j++) {
664 arrayRe[0][j]=(diagRe[0]*m.getRealElement(0,j) - diagIm[0]*m.getImagElement(0,j)) + (udiagRe[0]*m.getRealElement(1,j) - udiagIm[0]*m.getImagElement(1,j));
665 arrayIm[0][j]=(diagIm[0]*m.getRealElement(0,j) + diagRe[0]*m.getImagElement(0,j)) + (udiagIm[0]*m.getRealElement(1,j) + udiagRe[0]*m.getImagElement(1,j));
666 for(int i=1;i<lastRow;i++) {
667 arrayRe[i][j]=(ldiagRe[i]*m.getRealElement(i-1,j) - ldiagIm[i]*m.getImagElement(i-1,j)) + (diagRe[i]*m.getRealElement(i,j) - diagIm[i]*m.getImagElement(i,j)) + (udiagRe[i]*m.getRealElement(i+1,j) - udiagIm[i]*m.getImagElement(i+1,j));
668 arrayIm[i][j]=(ldiagIm[i]*m.getRealElement(i-1,j) + ldiagRe[i]*m.getImagElement(i-1,j)) + (diagIm[i]*m.getRealElement(i,j) + diagRe[i]*m.getImagElement(i,j)) + (udiagIm[i]*m.getRealElement(i+1,j) + udiagRe[i]*m.getImagElement(i+1,j));
669 }
670 arrayRe[lastRow][j]=(ldiagRe[lastRow]*m.getRealElement(lastRow-1,j) - ldiagIm[lastRow]*m.getImagElement(lastRow-1,j)) + (diagRe[lastRow]*m.getRealElement(lastRow,j) - diagIm[lastRow]*m.getImagElement(lastRow,j));
671 arrayIm[lastRow][j]=(ldiagIm[lastRow]*m.getRealElement(lastRow-1,j) + ldiagRe[lastRow]*m.getImagElement(lastRow-1,j)) + (diagIm[lastRow]*m.getRealElement(lastRow,j) + diagRe[lastRow]*m.getImagElement(lastRow,j));
672 }
673 return new ComplexSquareMatrix(arrayRe,arrayIm);
674 } else {
675 throw new MatrixDimensionException("Incompatible matrices.");
676 }
677 }
678 public ComplexSquareMatrix multiply(final ComplexSquareMatrix m) {
679 if(numCols==m.numRows) {
680 final double arrayRe[][]=new double[numRows][m.numCols];
681 final double arrayIm[][]=new double[numRows][m.numCols];
682 final int lastRow = numRows-1;
683 for(int j=0; j<numRows; j++) {
684 arrayRe[0][j]=(diagRe[0]*m.matrixRe[0][j] - diagIm[0]*m.matrixIm[0][j]) + (udiagRe[0]*m.matrixRe[1][j] - udiagIm[0]*m.matrixIm[1][j]);
685 arrayIm[0][j]=(diagIm[0]*m.matrixRe[0][j] + diagRe[0]*m.matrixIm[0][j]) + (udiagIm[0]*m.matrixRe[1][j] + udiagRe[0]*m.matrixIm[1][j]);
686 for(int i=1;i<lastRow;i++) {
687 arrayRe[i][j]=(ldiagRe[i]*m.matrixRe[i-1][j] - ldiagIm[i]*m.matrixIm[i-1][j]) + (diagRe[i]*m.matrixRe[i][j] - diagIm[i]*m.matrixIm[i][j]) + (udiagRe[i]*m.matrixRe[i+1][j] - udiagIm[i]*m.matrixIm[i+1][j]);
688 arrayIm[i][j]=(ldiagIm[i]*m.matrixRe[i-1][j] + ldiagRe[i]*m.matrixIm[i-1][j]) + (diagIm[i]*m.matrixRe[i][j] + diagRe[i]*m.matrixIm[i][j]) + (udiagIm[i]*m.matrixRe[i+1][j] + udiagRe[i]*m.matrixIm[i+1][j]);
689 }
690 arrayRe[lastRow][j]=(ldiagRe[lastRow]*m.matrixRe[lastRow-1][j] - ldiagIm[lastRow]*m.matrixIm[lastRow-1][j]) + (diagRe[lastRow]*m.matrixRe[lastRow][j] - diagIm[lastRow]*m.matrixIm[lastRow][j]);
691 arrayIm[lastRow][j]=(ldiagIm[lastRow]*m.matrixRe[lastRow-1][j] + ldiagRe[lastRow]*m.matrixIm[lastRow-1][j]) + (diagIm[lastRow]*m.matrixRe[lastRow][j] + diagRe[lastRow]*m.matrixIm[lastRow][j]);
692 }
693 return new ComplexSquareMatrix(arrayRe,arrayIm);
694 } else {
695 throw new MatrixDimensionException("Incompatible matrices.");
696 }
697 }
698
703 public ComplexSquareMatrix multiply(final ComplexTridiagonalMatrix m) {
704 int mRow=numRows;
705 if(numCols==m.numRows) {
706 final double arrayRe[][]=new double[mRow][mRow];
707 final double arrayIm[][]=new double[mRow][mRow];
708 arrayRe[0][0]=(diagRe[0]*m.diagRe[0] - diagIm[0]*m.diagIm[0])+(udiagRe[0]*m.ldiagRe[1] - udiagIm[0]*m.ldiagIm[1]);
709 arrayIm[0][0]=(diagIm[0]*m.diagRe[0] + diagRe[0]*m.diagIm[0])+(udiagIm[0]*m.ldiagRe[1] + udiagRe[0]*m.ldiagIm[1]);
710 arrayRe[0][1]=(diagRe[0]*m.udiagRe[0] - diagIm[0]*m.udiagIm[0])+(udiagRe[0]*m.diagRe[1] - udiagIm[0]*m.diagIm[1]);
711 arrayIm[0][1]=(diagIm[0]*m.udiagRe[0] + diagRe[0]*m.udiagIm[0])+(udiagIm[0]*m.diagRe[1] + udiagRe[0]*m.diagIm[1]);
712 arrayRe[0][2]=(udiagRe[0]*m.udiagRe[1] - udiagIm[0]*m.udiagIm[1]);
713 arrayIm[0][2]=(udiagIm[0]*m.udiagRe[1] + udiagRe[0]*m.udiagIm[1]);
714 if(mRow>3) {
715 arrayRe[1][0]=(ldiagRe[1]*m.diagRe[0] - ldiagIm[1]*m.diagIm[0])+(diagRe[1]*m.ldiagRe[1] - diagIm[1]*m.ldiagIm[1]);
716 arrayIm[1][0]=(ldiagIm[1]*m.diagRe[0] + ldiagRe[1]*m.diagIm[0])+(diagIm[1]*m.ldiagRe[1] + diagRe[1]*m.ldiagIm[1]);
717 arrayRe[1][1]=(ldiagRe[1]*m.udiagRe[0] - ldiagIm[1]*m.udiagIm[0])+(diagRe[1]*m.diagRe[1] - diagIm[1]*m.diagIm[1])+(udiagRe[1]*m.ldiagRe[2] - udiagIm[1]*m.ldiagIm[2]);
718 arrayIm[1][1]=(ldiagIm[1]*m.udiagRe[0] + ldiagRe[1]*m.udiagIm[0])+(diagIm[1]*m.diagRe[1] + diagRe[1]*m.diagIm[1])+(udiagIm[1]*m.ldiagRe[2] + udiagRe[1]*m.ldiagIm[2]);
719 arrayRe[1][2]=(diagRe[1]*m.udiagRe[1] - diagIm[1]*m.udiagIm[1])+(udiagRe[1]*m.diagRe[2] - udiagIm[1]*m.diagIm[2]);
720 arrayIm[1][2]=(diagIm[1]*m.udiagRe[1] + diagRe[1]*m.udiagIm[1])+(udiagIm[1]*m.diagRe[2] + udiagRe[1]*m.diagIm[2]);
721 arrayRe[1][3]=(udiagRe[1]*m.udiagRe[2] - udiagIm[1]*m.udiagIm[2]);
722 arrayIm[1][3]=(udiagIm[1]*m.udiagRe[2] + udiagRe[1]*m.udiagIm[2]);
723 }
724 if(mRow==3) {
725 arrayRe[1][0]=(ldiagRe[1]*m.diagRe[0] - ldiagIm[1]*m.diagIm[0])+(diagRe[1]*m.ldiagRe[1] - diagIm[1]*m.ldiagIm[1]);
726 arrayIm[1][0]=(ldiagIm[1]*m.diagRe[0] + ldiagRe[1]*m.diagIm[0])+(diagIm[1]*m.ldiagRe[1] + diagRe[1]*m.ldiagIm[1]);
727 arrayRe[1][1]=(ldiagRe[1]*m.udiagRe[0] - ldiagIm[1]*m.udiagIm[0])+(diagRe[1]*m.diagRe[1] - diagIm[1]*m.diagIm[1])+(udiagRe[1]*m.ldiagRe[2] - udiagIm[1]*m.ldiagIm[2]);
728 arrayIm[1][1]=(ldiagIm[1]*m.udiagRe[0] + ldiagRe[1]*m.udiagIm[0])+(diagIm[1]*m.diagRe[1] + diagRe[1]*m.diagIm[1])+(udiagIm[1]*m.ldiagRe[2] + udiagRe[1]*m.ldiagIm[2]);
729 arrayRe[1][2]=(diagRe[1]*m.udiagRe[1] - diagIm[1]*m.udiagIm[1])+(udiagRe[1]*m.diagRe[2] - udiagIm[1]*m.diagIm[2]);
730 arrayIm[1][2]=(diagIm[1]*m.udiagRe[1] + diagRe[1]*m.udiagIm[1])+(udiagIm[1]*m.diagRe[2] + udiagRe[1]*m.diagIm[2]);
731 } else if(mRow>4) {
732 for(int i=2;i<mRow-2;i++) {
733 arrayRe[i][i-2]=(ldiagRe[i]*m.ldiagRe[i-1] - ldiagIm[i]*m.ldiagIm[i-1]);
734 arrayIm[i][i-2]=(ldiagIm[i]*m.ldiagRe[i-1] + ldiagRe[i]*m.ldiagIm[i-1]);
735 arrayRe[i][i-1]=(ldiagRe[i]*m.diagRe[i-1] - ldiagIm[i]*m.diagIm[i-1])+(diagRe[i]*m.ldiagRe[i] - diagIm[i]*m.ldiagIm[i]);
736 arrayIm[i][i-1]=(ldiagIm[i]*m.diagRe[i-1] + ldiagRe[i]*m.diagIm[i-1])+(diagIm[i]*m.ldiagRe[i] + diagRe[i]*m.ldiagIm[i]);
737 arrayRe[i][i]=(ldiagRe[i]*m.udiagRe[i-1] - ldiagIm[i]*m.udiagIm[i-1])+(diagRe[i]*m.diagRe[i] - diagIm[i]*m.diagIm[i])+(udiagRe[i]*m.ldiagRe[i+1] - udiagIm[i]*m.ldiagIm[i+1]);
738 arrayIm[i][i]=(ldiagIm[i]*m.udiagRe[i-1] + ldiagRe[i]*m.udiagIm[i-1])+(diagIm[i]*m.diagRe[i] + diagRe[i]*m.diagIm[i])+(udiagIm[i]*m.ldiagRe[i+1] + udiagRe[i]*m.ldiagIm[i+1]);
739 arrayRe[i][i+1]=(diagRe[i]*m.udiagRe[i] - diagIm[i]*m.udiagIm[i])+(udiagRe[i]*m.diagRe[i+1] - udiagIm[i]*m.diagIm[i+1]);
740 arrayIm[i][i+1]=(diagIm[i]*m.udiagRe[i] + diagRe[i]*m.udiagIm[i])+(udiagIm[i]*m.diagRe[i+1] + udiagRe[i]*m.diagIm[i+1]);
741 arrayRe[i][i+2]=(udiagRe[i]*m.udiagRe[i+1] - udiagIm[i]*m.udiagIm[i+1]);
742 arrayIm[i][i+2]=(udiagIm[i]*m.udiagRe[i+1] + udiagRe[i]*m.udiagIm[i+1]);
743 }
744 }
745 if(mRow>3) {
746 arrayRe[mRow-2][mRow-4]=(ldiagRe[mRow-2]*m.ldiagRe[mRow-3] - ldiagIm[mRow-2]*m.ldiagIm[mRow-3]);
747 arrayIm[mRow-2][mRow-4]=(ldiagIm[mRow-2]*m.ldiagRe[mRow-3] + ldiagRe[mRow-2]*m.ldiagIm[mRow-3]);
748 arrayRe[mRow-2][mRow-3]=(ldiagRe[mRow-2]*m.diagRe[mRow-3] - ldiagIm[mRow-2]*m.diagIm[mRow-3])+(diagRe[mRow-2]*m.ldiagRe[mRow-2] - diagIm[mRow-2]*m.ldiagIm[mRow-2]);
749 arrayIm[mRow-2][mRow-3]=(ldiagIm[mRow-2]*m.diagRe[mRow-3] + ldiagRe[mRow-2]*m.diagIm[mRow-3])+(diagIm[mRow-2]*m.ldiagRe[mRow-2] + diagRe[mRow-2]*m.ldiagIm[mRow-2]);
750 arrayRe[mRow-2][mRow-2]=(ldiagRe[mRow-2]*m.udiagRe[mRow-3] - ldiagIm[mRow-2]*m.udiagIm[mRow-3])+(diagRe[mRow-2]*m.diagRe[mRow-2] - diagIm[mRow-2]*m.diagIm[mRow-2])+(udiagRe[mRow-2]*m.ldiagRe[mRow-1] - udiagIm[mRow-2]*m.ldiagIm[mRow-1]);
751 arrayIm[mRow-2][mRow-2]=(ldiagIm[mRow-2]*m.udiagRe[mRow-3] + ldiagRe[mRow-2]*m.udiagIm[mRow-3])+(diagIm[mRow-2]*m.diagRe[mRow-2] + diagRe[mRow-2]*m.diagIm[mRow-2])+(udiagIm[mRow-2]*m.ldiagRe[mRow-1] + udiagRe[mRow-2]*m.ldiagIm[mRow-1]);
752 arrayRe[mRow-2][mRow-1]=(diagRe[mRow-2]*m.udiagRe[mRow-2] - diagIm[mRow-2]*m.udiagIm[mRow-2])+(udiagRe[mRow-2]*m.diagRe[mRow-1] - udiagIm[mRow-2]*m.diagIm[mRow-1]);
753 arrayIm[mRow-2][mRow-1]=(diagIm[mRow-2]*m.udiagRe[mRow-2] + diagRe[mRow-2]*m.udiagIm[mRow-2])+(udiagIm[mRow-2]*m.diagRe[mRow-1] + udiagRe[mRow-2]*m.diagIm[mRow-1]);
754 }
755 mRow--;
756 arrayRe[mRow][mRow-2]=(ldiagRe[mRow]*m.ldiagRe[mRow-1] - ldiagIm[mRow]*m.ldiagIm[mRow-1]);
757 arrayIm[mRow][mRow-2]=(ldiagIm[mRow]*m.ldiagRe[mRow-1] + ldiagRe[mRow]*m.ldiagIm[mRow-1]);
758 arrayRe[mRow][mRow-1]=(ldiagRe[mRow]*m.diagRe[mRow-1] - ldiagIm[mRow]*m.diagIm[mRow-1])+(diagRe[mRow]*m.ldiagRe[mRow] - diagIm[mRow]*m.ldiagIm[mRow]);
759 arrayIm[mRow][mRow-1]=(ldiagIm[mRow]*m.diagRe[mRow-1] + ldiagRe[mRow]*m.diagIm[mRow-1])+(diagIm[mRow]*m.ldiagRe[mRow] + diagRe[mRow]*m.ldiagIm[mRow]);
760 arrayRe[mRow][mRow]=(ldiagRe[mRow]*m.udiagRe[mRow-1] - ldiagIm[mRow]*m.udiagIm[mRow-1])+(diagRe[mRow]*m.diagRe[mRow] - diagIm[mRow]*m.diagIm[mRow]);
761 arrayIm[mRow][mRow]=(ldiagIm[mRow]*m.udiagRe[mRow-1] + ldiagRe[mRow]*m.udiagIm[mRow-1])+(diagIm[mRow]*m.diagRe[mRow] + diagRe[mRow]*m.diagIm[mRow]);
762 return new ComplexSquareMatrix(arrayRe,arrayIm);
763 } else
764 throw new MatrixDimensionException("Incompatible matrices.");
765 }
766 private ComplexSquareMatrix multiplyTridiagonal(final AbstractComplexSquareMatrix m) {
767 int mRow=numRows;
768 if(numCols==m.rows()) {
769 final double arrayRe[][]=new double[mRow][mRow];
770 final double arrayIm[][]=new double[mRow][mRow];
771 Complex elem1,elem2,elem3;
772 elem1=m.getElement(0,0);elem2=m.getElement(1,0);
773 arrayRe[0][0]=(diagRe[0]*elem1.real() - diagIm[0]*elem1.imag())+(udiagRe[0]*elem2.real() - udiagIm[0]*elem2.imag());
774 arrayIm[0][0]=(diagIm[0]*elem1.real() + diagRe[0]*elem1.imag())+(udiagIm[0]*elem2.real() + udiagRe[0]*elem2.imag());
775 elem1=m.getElement(0,1);elem2=m.getElement(1,1);
776 arrayRe[0][1]=(diagRe[0]*elem1.real() - diagIm[0]*elem1.imag())+(udiagRe[0]*elem2.real() - udiagIm[0]*elem2.imag());
777 arrayIm[0][1]=(diagIm[0]*elem1.real() + diagRe[0]*elem1.imag())+(udiagIm[0]*elem2.real() + udiagRe[0]*elem2.imag());
778 elem1=m.getElement(1,2);
779 arrayRe[0][2]=(udiagRe[0]*elem1.real() - udiagIm[0]*elem1.imag());
780 arrayIm[0][2]=(udiagIm[0]*elem1.real() + udiagRe[0]*elem1.imag());
781 if(mRow>3) {
782 elem1=m.getElement(0,0);elem2=m.getElement(1,0);
783 arrayRe[1][0]=(ldiagRe[1]*elem1.real() - ldiagIm[1]*elem1.imag())+(diagRe[1]*elem2.real() - diagIm[1]*elem2.imag());
784 arrayIm[1][0]=(ldiagIm[1]*elem1.real() + ldiagRe[1]*elem1.imag())+(diagIm[1]*elem2.real() + diagRe[1]*elem2.imag());
785 elem1=m.getElement(0,1);elem2=m.getElement(1,1);elem3=m.getElement(2,1);
786 arrayRe[1][1]=(ldiagRe[1]*elem1.real() - ldiagIm[1]*elem1.imag())+(diagRe[1]*elem2.real() - diagIm[1]*elem2.imag())+(udiagRe[1]*elem3.real() - udiagIm[1]*elem3.imag());
787 arrayIm[1][1]=(ldiagIm[1]*elem1.real() + ldiagRe[1]*elem1.imag())+(diagIm[1]*elem2.real() + diagRe[1]*elem2.imag())+(udiagIm[1]*elem3.real() + udiagRe[1]*elem3.imag());
788 elem1=m.getElement(1,2);elem2=m.getElement(2,2);
789 arrayRe[1][2]=(diagRe[1]*elem1.real() - diagIm[1]*elem1.imag())+(udiagRe[1]*elem2.real() - udiagIm[1]*elem2.imag());
790 arrayIm[1][2]=(diagIm[1]*elem1.real() + diagRe[1]*elem1.imag())+(udiagIm[1]*elem2.real() + udiagRe[1]*elem2.imag());
791 elem1=m.getElement(2,3);
792 arrayRe[1][3]=(udiagRe[1]*elem1.real() - udiagIm[1]*elem1.imag());
793 arrayIm[1][3]=(udiagIm[1]*elem1.real() + udiagRe[1]*elem1.imag());
794 }
795 if(mRow==3) {
796 elem1=m.getElement(0,0);elem2=m.getElement(1,0);
797 arrayRe[1][0]=(ldiagRe[1]*elem1.real() - ldiagIm[1]*elem1.imag())+(diagRe[1]*elem2.real() - diagIm[1]*elem2.imag());
798 arrayIm[1][0]=(ldiagIm[1]*elem1.real() + ldiagRe[1]*elem1.imag())+(diagIm[1]*elem2.real() + diagRe[1]*elem2.imag());
799 elem1=m.getElement(0,1);elem2=m.getElement(1,1);elem3=m.getElement(2,1);
800 arrayRe[1][1]=(ldiagRe[1]*elem1.real() - ldiagIm[1]*elem1.imag())+(diagRe[1]*elem2.real() - diagIm[1]*elem2.imag())+(udiagRe[1]*elem3.real() - udiagIm[1]*elem3.imag());
801 arrayIm[1][1]=(ldiagIm[1]*elem1.real() + ldiagRe[1]*elem1.imag())+(diagIm[1]*elem2.real() + diagRe[1]*elem2.imag())+(udiagIm[1]*elem3.real() + udiagRe[1]*elem3.imag());
802 elem1=m.getElement(1,2);elem2=m.getElement(2,2);
803 arrayRe[1][2]=(diagRe[1]*elem1.real() - diagIm[1]*elem1.imag())+(udiagRe[1]*elem2.real() - udiagIm[1]*elem2.imag());
804 arrayIm[1][2]=(diagIm[1]*elem1.real() + diagRe[1]*elem1.imag())+(udiagIm[1]*elem2.real() + udiagRe[1]*elem2.imag());
805 } else if(mRow>4) {
806 for(int i=2;i<mRow-2;i++) {
807 elem1=m.getElement(i-1,i-2);
808 arrayRe[i][i-2]=(ldiagRe[i]*elem1.real() - ldiagIm[i]*elem1.imag());
809 arrayIm[i][i-2]=(ldiagIm[i]*elem1.real() + ldiagRe[i]*elem1.imag());
810 elem1=m.getElement(i-1,i-1);elem2=m.getElement(i,i-1);
811 arrayRe[i][i-1]=(ldiagRe[i]*elem1.real() - ldiagIm[i]*elem1.imag())+(diagRe[i]*elem2.real() - diagIm[i]*elem2.imag());
812 arrayIm[i][i-1]=(ldiagIm[i]*elem1.real() + ldiagRe[i]*elem1.imag())+(diagIm[i]*elem2.real() + diagRe[i]*elem2.imag());
813 elem1=m.getElement(i-1,i);elem2=m.getElement(i,i);elem3=m.getElement(i+1,i);
814 arrayRe[i][i]=(ldiagRe[i]*elem1.real() - ldiagIm[i]*elem1.imag())+(diagRe[i]*elem2.real() - diagIm[i]*elem2.imag())+(udiagRe[i]*elem3.real() - udiagIm[i]*elem3.imag());
815 arrayIm[i][i]=(ldiagIm[i]*elem1.real() + ldiagRe[i]*elem1.imag())+(diagIm[i]*elem2.real() + diagRe[i]*elem2.imag())+(udiagIm[i]*elem3.real() + udiagRe[i]*elem3.imag());
816 elem1=m.getElement(i,i+1);elem2=m.getElement(i+1,i+1);
817 arrayRe[i][i+1]=(diagRe[i]*elem1.real() - diagIm[i]*elem1.imag())+(udiagRe[i]*elem2.real() - udiagIm[i]*elem2.imag());
818 arrayIm[i][i+1]=(diagIm[i]*elem1.real() + diagRe[i]*elem1.imag())+(udiagIm[i]*elem2.real() + udiagRe[i]*elem2.imag());
819 elem1=m.getElement(i+1,i+2);
820 arrayRe[i][i+2]=(udiagRe[i]*elem1.real() - udiagIm[i]*elem1.imag());
821 arrayIm[i][i+2]=(udiagIm[i]*elem1.real() + udiagRe[i]*elem1.imag());
822 }
823 }
824 if(mRow>3) {
825 elem1=m.getElement(mRow-3,mRow-4);
826 arrayRe[mRow-2][mRow-4]=(ldiagRe[mRow-2]*elem1.real() - ldiagIm[mRow-2]*elem1.imag());
827 arrayIm[mRow-2][mRow-4]=(ldiagIm[mRow-2]*elem1.real() + ldiagRe[mRow-2]*elem1.imag());
828 elem1=m.getElement(mRow-3,mRow-3);elem2=m.getElement(mRow-2,mRow-3);
829 arrayRe[mRow-2][mRow-3]=(ldiagRe[mRow-2]*elem1.real() - ldiagIm[mRow-2]*elem1.imag())+(diagRe[mRow-2]*elem2.real() - diagIm[mRow-2]*elem2.imag());
830 arrayIm[mRow-2][mRow-3]=(ldiagIm[mRow-2]*elem1.real() + ldiagRe[mRow-2]*elem1.imag())+(diagIm[mRow-2]*elem2.real() + diagRe[mRow-2]*elem2.imag());
831 elem1=m.getElement(mRow-3,mRow-2);elem2=m.getElement(mRow-2,mRow-2);elem3=m.getElement(mRow-1,mRow-2);
832 arrayRe[mRow-2][mRow-2]=(ldiagRe[mRow-2]*elem1.real() - ldiagIm[mRow-2]*elem1.imag())+(diagRe[mRow-2]*elem2.real() - diagIm[mRow-2]*elem2.imag())+(udiagRe[mRow-2]*elem3.real() - udiagIm[mRow-2]*elem3.imag());
833 arrayIm[mRow-2][mRow-2]=(ldiagIm[mRow-2]*elem1.real() + ldiagRe[mRow-2]*elem1.imag())+(diagIm[mRow-2]*elem2.real() + diagRe[mRow-2]*elem2.imag())+(udiagIm[mRow-2]*elem3.real() + udiagRe[mRow-2]*elem3.imag());
834 elem1=m.getElement(mRow-2,mRow-1);elem2=m.getElement(mRow-1,mRow-1);
835 arrayRe[mRow-2][mRow-1]=(diagRe[mRow-2]*elem1.real() - diagIm[mRow-2]*elem1.imag())+(udiagRe[mRow-2]*elem2.real() - udiagIm[mRow-2]*elem2.imag());
836 arrayIm[mRow-2][mRow-1]=(diagIm[mRow-2]*elem1.real() + diagRe[mRow-2]*elem1.imag())+(udiagIm[mRow-2]*elem2.real() + udiagRe[mRow-2]*elem2.imag());
837 }
838 mRow--;
839 elem1=m.getElement(mRow-1,mRow-2);
840 arrayRe[mRow][mRow-2]=(ldiagRe[mRow]*elem1.real() - ldiagIm[mRow]*elem1.imag());
841 arrayIm[mRow][mRow-2]=(ldiagIm[mRow]*elem1.real() + ldiagRe[mRow]*elem1.imag());
842 elem1=m.getElement(mRow-1,mRow-1);elem2=m.getElement(mRow,mRow-1);
843 arrayRe[mRow][mRow-1]=(ldiagRe[mRow]*elem1.real() - ldiagIm[mRow]*elem1.imag())+(diagRe[mRow]*elem2.real() - diagIm[mRow]*elem2.imag());
844 arrayIm[mRow][mRow-1]=(ldiagIm[mRow]*elem1.real() + ldiagRe[mRow]*elem1.imag())+(diagIm[mRow]*elem2.real() + diagRe[mRow]*elem2.imag());
845 elem1=m.getElement(mRow-1,mRow);elem2=m.getElement(mRow,mRow);
846 arrayRe[mRow][mRow]=(ldiagRe[mRow]*elem1.real() - ldiagIm[mRow]*elem1.imag())+(diagRe[mRow]*elem2.real() - diagIm[mRow]*elem2.imag());
847 arrayIm[mRow][mRow]=(ldiagIm[mRow]*elem1.real() + ldiagRe[mRow]*elem1.imag())+(diagIm[mRow]*elem2.real() + diagRe[mRow]*elem2.imag());
848 return new ComplexSquareMatrix(arrayRe,arrayIm);
849 } else {
850 throw new MatrixDimensionException("Incompatible matrices.");
851 }
852 }
853
854
856
860 public AbstractComplexMatrix hermitianAdjoint() {
861 int mRow=numRows;
862 final ComplexTridiagonalMatrix ans=new ComplexTridiagonalMatrix(mRow);
863 System.arraycopy(ldiagRe,1,ans.udiagRe,0,ldiagRe.length-1);
864 System.arraycopy(diagRe,0,ans.diagRe,0,diagRe.length);
865 System.arraycopy(udiagRe,0,ans.ldiagRe,1,udiagRe.length-1);
866 ans.diagIm[0]=-diagIm[0];
867 ans.ldiagIm[1]=-udiagIm[0];
868 mRow--;
869 for(int i=1;i<mRow;i++) {
870 ans.udiagIm[i-1]=-ldiagIm[i];
871 ans.diagIm[i]=-diagIm[i];
872 ans.ldiagIm[i+1]=-udiagIm[i];
873 }
874 ans.udiagIm[mRow-1]=-ldiagIm[mRow];
875 ans.diagIm[mRow]=-diagIm[mRow];
876 return ans;
877 }
878
879
881
885 public AbstractComplexMatrix conjugate() {
886 int mRow=numRows;
887 final ComplexTridiagonalMatrix ans=new ComplexTridiagonalMatrix(mRow);
888 System.arraycopy(ldiagRe,1,ans.ldiagRe,0,ldiagRe.length-1);
889 System.arraycopy(diagRe,0,ans.diagRe,0,diagRe.length);
890 System.arraycopy(udiagRe,0,ans.udiagRe,1,udiagRe.length-1);
891 ans.diagIm[0]=-diagIm[0];
892 ans.udiagIm[0]=-udiagIm[0];
893 mRow--;
894 for(int i=1;i<mRow;i++) {
895 ans.ldiagIm[i]=-ldiagIm[i];
896 ans.diagIm[i]=-diagIm[i];
897 ans.udiagIm[i]=-udiagIm[i];
898 }
899 ans.ldiagIm[mRow]=-ldiagIm[mRow];
900 ans.diagIm[mRow]=-diagIm[mRow];
901 return ans;
902 }
903
904
906
910 public Matrix transpose() {
911 final ComplexTridiagonalMatrix ans=new ComplexTridiagonalMatrix(numRows);
912 System.arraycopy(ldiagRe,1,ans.udiagRe,0,ldiagRe.length-1);
913 System.arraycopy(ldiagIm,1,ans.udiagIm,0,ldiagIm.length-1);
914 System.arraycopy(diagRe,0,ans.diagRe,0,diagRe.length);
915 System.arraycopy(diagIm,0,ans.diagIm,0,diagIm.length);
916 System.arraycopy(udiagRe,0,ans.ldiagRe,1,udiagRe.length-1);
917 System.arraycopy(udiagIm,0,ans.ldiagIm,1,udiagIm.length-1);
918 return ans;
919 }
920
921
923
928 public AbstractComplexMatrix mapElements(final ComplexMapping f) {
929 Complex zeroValue = f.map(Complex.ZERO);
930 if(zeroValue.mod() <= JSci.GlobalSettings.ZERO_TOL)
931 return tridiagonalMap(f);
932 else
933 return generalMap(f, zeroValue);
934 }
935 private AbstractComplexMatrix tridiagonalMap(ComplexMapping f) {
936 int mRow=numRows;
937 final ComplexTridiagonalMatrix ans=new ComplexTridiagonalMatrix(mRow);
938 ans.setElement(0,0,f.map(diagRe[0],diagIm[0]));
939 ans.setElement(0,1,f.map(udiagRe[0],udiagIm[0]));
940 mRow--;
941 for(int i=1;i<mRow;i++) {
942 ans.setElement(i,i-1,f.map(ldiagRe[i],ldiagIm[i]));
943 ans.setElement(i,i,f.map(diagRe[i],diagIm[i]));
944 ans.setElement(i,i+1,f.map(udiagRe[i],udiagIm[i]));
945 }
946 ans.setElement(mRow,mRow-1,f.map(ldiagRe[mRow],ldiagIm[mRow]));
947 ans.setElement(mRow,mRow,f.map(diagRe[mRow],diagIm[mRow]));
948 return ans;
949 }
950 private AbstractComplexMatrix generalMap(ComplexMapping f, Complex zeroValue) {
951 final double arrayRe[][]=new double[numRows][numRows];
952 final double arrayIm[][]=new double[numRows][numRows];
953 for(int i=0; i<numRows; i++) {
954 for(int j=0; j<numRows; j++) {
955 arrayRe[i][j] = zeroValue.real();
956 arrayIm[i][j] = zeroValue.imag();
957 }
958 }
959 int mRow=numRows;
960 Complex z = f.map(diagRe[0], diagIm[0]);
961 arrayRe[0][0]=z.real(); arrayIm[0][0]=z.imag();
962 z = f.map(udiagRe[0], udiagIm[0]);
963 arrayRe[0][1]=z.real(); arrayIm[0][1]=z.imag();
964 mRow--;
965 for(int i=1;i<mRow;i++) {
966 z = f.map(ldiagRe[i], ldiagIm[i]);
967 arrayRe[i][i-1]=z.real(); arrayIm[i][i-1]=z.imag();
968 z = f.map(diagRe[i], diagIm[i]);
969 arrayRe[i][i]=z.real(); arrayIm[i][i]=z.imag();
970 z = f.map(udiagRe[i], udiagIm[i]);
971 arrayRe[i][i+1]=z.real(); arrayIm[i][i+1]=z.imag();
972 }
973 z = f.map(ldiagRe[mRow], ldiagIm[mRow]);
974 arrayRe[mRow][mRow-1]=z.real(); arrayIm[mRow][mRow-1]=z.imag();
975 z = f.map(diagRe[mRow], diagIm[mRow]);
976 arrayRe[mRow][mRow]=z.real(); arrayIm[mRow][mRow]=z.imag();
977 return new ComplexSquareMatrix(arrayRe, arrayIm);
978 }
979 }
980 |
Popular Tags |
Java API By Example, From Geeks To Geeks.