0001: /* AUTO-GENERATED */
0002: package JSci.maths.matrices;
0003:
0004: import JSci.maths.ArrayMath;
0005: import JSci.maths.Complex;
0006: import JSci.maths.ComplexMapping;
0007: import JSci.maths.LinearMath;
0008: import JSci.maths.DimensionException;
0009: import JSci.maths.MaximumIterationsExceededException;
0010: import JSci.maths.vectors.AbstractComplexVector;
0011: import JSci.maths.vectors.ComplexVector;
0012:
0013: /**
0014: * The ComplexTridiagonalMatrix class provides an object for encapsulating tridiagonal matrices containing complex numbers.
0015: * Uses compressed diagonal storage.
0016: * @version 2.2
0017: * @author Mark Hale
0018: */
0019: public class ComplexTridiagonalMatrix extends
0020: AbstractComplexSquareMatrix implements TridiagonalMatrix {
0021: /**
0022: * Tridiagonal data.
0023: */
0024: protected final double ldiagRe[], ldiagIm[];
0025: protected final double diagRe[], diagIm[];
0026: protected final double udiagRe[], udiagIm[];
0027:
0028: /**
0029: * Constructs an empty matrix.
0030: * @param size the number of rows/columns
0031: */
0032: public ComplexTridiagonalMatrix(final int size) {
0033: super (size);
0034: ldiagRe = new double[size];
0035: ldiagIm = new double[size];
0036: diagRe = new double[size];
0037: diagIm = new double[size];
0038: udiagRe = new double[size];
0039: udiagIm = new double[size];
0040: }
0041:
0042: /**
0043: * Constructs a matrix from an array.
0044: * Any non-tridiagonal elements in the array are ignored.
0045: * @param array an assigned value
0046: * @exception MatrixDimensionException If the array is not square.
0047: */
0048: public ComplexTridiagonalMatrix(final Complex array[][]) {
0049: this (array.length);
0050: if (!ArrayMath.isSquare(array))
0051: throw new MatrixDimensionException("Array is not square.");
0052: diagRe[0] = array[0][0].real();
0053: diagIm[0] = array[0][0].imag();
0054: udiagRe[0] = array[0][1].real();
0055: udiagIm[0] = array[0][1].imag();
0056: int i = 1;
0057: for (; i < array.length - 1; i++) {
0058: ldiagRe[i] = array[i][i - 1].real();
0059: ldiagIm[i] = array[i][i - 1].imag();
0060: diagRe[i] = array[i][i].real();
0061: diagIm[i] = array[i][i].imag();
0062: udiagRe[i] = array[i][i + 1].real();
0063: udiagIm[i] = array[i][i + 1].imag();
0064: }
0065: ldiagRe[i] = array[i][i - 1].real();
0066: ldiagIm[i] = array[i][i - 1].imag();
0067: diagRe[i] = array[i][i].real();
0068: diagIm[i] = array[i][i].imag();
0069: }
0070:
0071: /**
0072: * Compares two complex matrices for equality.
0073: * @param m a complex matrix
0074: */
0075: public boolean equals(AbstractComplexMatrix m, double tol) {
0076: if (m instanceof TridiagonalMatrix) {
0077: if (numRows != m.rows() || numCols != m.columns())
0078: return false;
0079: double sumSqr = 0;
0080: double ldeltaRe, ldeltaIm;
0081: double deltaRe = diagRe[0] - m.getRealElement(0, 0);
0082: double deltaIm = diagIm[0] - m.getImagElement(0, 0);
0083: double udeltaRe = udiagRe[0] - m.getRealElement(0, 1);
0084: double udeltaIm = udiagIm[0] - m.getImagElement(0, 1);
0085: sumSqr += deltaRe * deltaRe + deltaIm * deltaIm + udeltaRe
0086: * udeltaRe + udeltaIm * udeltaIm;
0087: int i = 1;
0088: for (; i < numRows - 1; i++) {
0089: ldeltaRe = ldiagRe[i] - m.getRealElement(i, i - 1);
0090: ldeltaIm = ldiagIm[i] - m.getImagElement(i, i - 1);
0091: deltaRe = diagRe[i] - m.getRealElement(i, i);
0092: deltaIm = diagIm[i] - m.getImagElement(i, i);
0093: udeltaRe = udiagRe[i] - m.getRealElement(i, i + 1);
0094: udeltaIm = udiagIm[i] - m.getImagElement(i, i + 1);
0095: sumSqr += ldeltaRe * ldeltaRe + ldeltaIm * ldeltaIm
0096: + deltaRe * deltaRe + deltaIm * deltaIm
0097: + udeltaRe * udeltaRe + udeltaIm * udeltaIm;
0098: }
0099: ldeltaRe = ldiagRe[i] - m.getRealElement(i, i - 1);
0100: ldeltaIm = ldiagIm[i] - m.getImagElement(i, i - 1);
0101: deltaRe = diagRe[i] - m.getRealElement(i, i);
0102: deltaIm = diagIm[i] - m.getImagElement(i, i);
0103: sumSqr += ldeltaRe * ldeltaRe + ldeltaIm * ldeltaIm
0104: + deltaRe * deltaRe + deltaIm * deltaIm;
0105: return (sumSqr <= tol * tol);
0106: } else {
0107: return false;
0108: }
0109: }
0110:
0111: /**
0112: * Returns a string representing this matrix.
0113: */
0114: public String toString() {
0115: final StringBuffer buf = new StringBuffer(5 * rows()
0116: * columns());
0117: for (int i = 0; i < rows(); i++) {
0118: for (int j = 0; j < columns(); j++) {
0119: buf.append(getElement(i, j).toString());
0120: buf.append(' ');
0121: }
0122: buf.append('\n');
0123: }
0124: return buf.toString();
0125: }
0126:
0127: /**
0128: * Returns the real part of this complex matrix.
0129: * @return a double tridiagonal matrix
0130: */
0131: public AbstractDoubleMatrix real() {
0132: final DoubleTridiagonalMatrix m = new DoubleTridiagonalMatrix(
0133: numRows);
0134: m.diag[0] = diagRe[0];
0135: m.udiag[0] = udiagRe[0];
0136: int i = 1;
0137: for (; i < numRows - 1; i++) {
0138: m.ldiag[i] = ldiagRe[i];
0139: m.diag[i] = diagRe[i];
0140: m.udiag[i] = udiagRe[i];
0141: }
0142: m.ldiag[i] = ldiagRe[i];
0143: m.diag[i] = diagRe[i];
0144: return m;
0145: }
0146:
0147: /**
0148: * Returns the imaginary part of this complex matrix.
0149: * @return a double tridiagonal matrix
0150: */
0151: public AbstractDoubleMatrix imag() {
0152: final DoubleTridiagonalMatrix m = new DoubleTridiagonalMatrix(
0153: numRows);
0154: m.diag[0] = diagIm[0];
0155: m.udiag[0] = udiagIm[0];
0156: int i = 1;
0157: for (; i < numRows - 1; i++) {
0158: m.ldiag[i] = ldiagIm[i];
0159: m.diag[i] = diagIm[i];
0160: m.udiag[i] = udiagIm[i];
0161: }
0162: m.ldiag[i] = ldiagIm[i];
0163: m.diag[i] = diagIm[i];
0164: return m;
0165: }
0166:
0167: /**
0168: * Returns an element of the matrix.
0169: * @param i row index of the element
0170: * @param j column index of the element
0171: * @exception MatrixDimensionException If attempting to access an invalid element.
0172: */
0173: public Complex getElement(final int i, final int j) {
0174: if (i >= 0 && i < numRows && j >= 0 && j < numCols) {
0175: if (j == i - 1)
0176: return new Complex(ldiagRe[i], ldiagIm[i]);
0177: else if (j == i)
0178: return new Complex(diagRe[i], diagIm[i]);
0179: else if (j == i + 1)
0180: return new Complex(udiagRe[i], udiagIm[i]);
0181: else
0182: return Complex.ZERO;
0183: } else
0184: throw new MatrixDimensionException(getInvalidElementMsg(i,
0185: j));
0186: }
0187:
0188: public double getRealElement(final int i, final int j) {
0189: if (i >= 0 && i < numRows && j >= 0 && j < numCols) {
0190: if (j == i - 1)
0191: return ldiagRe[i];
0192: else if (j == i)
0193: return diagRe[i];
0194: else if (j == i + 1)
0195: return udiagRe[i];
0196: else
0197: return 0.0;
0198: } else
0199: throw new MatrixDimensionException(getInvalidElementMsg(i,
0200: j));
0201: }
0202:
0203: public double getImagElement(final int i, final int j) {
0204: if (i >= 0 && i < numRows && j >= 0 && j < numCols) {
0205: if (j == i - 1)
0206: return ldiagIm[i];
0207: else if (j == i)
0208: return diagIm[i];
0209: else if (j == i + 1)
0210: return udiagIm[i];
0211: else
0212: return 0.0;
0213: } else
0214: throw new MatrixDimensionException(getInvalidElementMsg(i,
0215: j));
0216: }
0217:
0218: /**
0219: * Sets the value of an element of the matrix.
0220: * Should only be used to initialise this matrix.
0221: * @param i row index of the element
0222: * @param j column index of the element
0223: * @param z a complex number
0224: * @exception MatrixDimensionException If attempting to access an invalid element.
0225: */
0226: public void setElement(final int i, final int j, final Complex z) {
0227: if (i >= 0 && i < numRows && j >= 0 && j < numCols) {
0228: if (j == i - 1) {
0229: ldiagRe[i] = z.real();
0230: ldiagIm[i] = z.imag();
0231: } else if (j == i) {
0232: diagRe[i] = z.real();
0233: diagIm[i] = z.imag();
0234: } else if (j == i + 1) {
0235: udiagRe[i] = z.real();
0236: udiagIm[i] = z.imag();
0237: } else {
0238: throw new MatrixDimensionException(
0239: getInvalidElementMsg(i, j));
0240: }
0241: } else
0242: throw new MatrixDimensionException(getInvalidElementMsg(i,
0243: j));
0244: }
0245:
0246: /**
0247: * Sets the value of an element of the matrix.
0248: * Should only be used to initialise this matrix.
0249: * @param i row index of the element
0250: * @param j column index of the element
0251: * @param x the real part of a complex number
0252: * @param y the imaginary part of a complex number
0253: * @exception MatrixDimensionException If attempting to access an invalid element.
0254: */
0255: public void setElement(final int i, final int j, final double x,
0256: final double y) {
0257: if (i >= 0 && i < numRows && j >= 0 && j < numCols) {
0258: if (j == i - 1) {
0259: ldiagRe[i] = x;
0260: ldiagIm[i] = y;
0261: } else if (j == i) {
0262: diagRe[i] = x;
0263: diagIm[i] = y;
0264: } else if (j == i + 1) {
0265: udiagRe[i] = x;
0266: udiagIm[i] = y;
0267: } else {
0268: throw new MatrixDimensionException(
0269: getInvalidElementMsg(i, j));
0270: }
0271: } else
0272: throw new MatrixDimensionException(getInvalidElementMsg(i,
0273: j));
0274: }
0275:
0276: /**
0277: * Returns the trace.
0278: */
0279: public Complex trace() {
0280: double trRe = diagRe[0];
0281: double trIm = diagIm[0];
0282: for (int i = 1; i < numRows; i++) {
0283: trRe += diagRe[i];
0284: trIm += diagIm[i];
0285: }
0286: return new Complex(trRe, trIm);
0287: }
0288:
0289: /**
0290: * Returns the l<sup><img border=0 alt="infinity" src="doc-files/infinity.gif"></sup>-norm.
0291: * @author Taber Smith
0292: */
0293: public double infNorm() {
0294: double result = Math.sqrt((diagRe[0] * diagRe[0] + diagIm[0]
0295: * diagIm[0]))
0296: + Math.sqrt((udiagRe[0] * udiagRe[0] + udiagIm[0]
0297: * udiagIm[0]));
0298: double tmpResult;
0299: int i = 1;
0300: for (; i < numRows - 1; i++) {
0301: tmpResult = Math.sqrt((ldiagRe[i] * ldiagRe[i] + ldiagIm[i]
0302: * ldiagIm[i]))
0303: + Math.sqrt((diagRe[i] * diagRe[i] + diagIm[i]
0304: * diagIm[i]))
0305: + Math.sqrt((udiagRe[i] * udiagRe[i] + udiagIm[i]
0306: * udiagIm[i]));
0307: if (tmpResult > result)
0308: result = tmpResult;
0309: }
0310: tmpResult = Math.sqrt((ldiagRe[i] * ldiagRe[i] + ldiagIm[i]
0311: * ldiagIm[i]))
0312: + Math.sqrt((diagRe[i] * diagRe[i] + diagIm[i]
0313: * diagIm[i]));
0314: if (tmpResult > result)
0315: result = tmpResult;
0316: return result;
0317: }
0318:
0319: /**
0320: * Returns the Frobenius (l<sup>2</sup>) norm.
0321: * @author Taber Smith
0322: */
0323: public double frobeniusNorm() {
0324: double result = diagRe[0] * diagRe[0] + diagIm[0] * diagIm[0]
0325: + udiagRe[0] * udiagRe[0] + udiagIm[0] * udiagIm[0];
0326: int i = 1;
0327: for (; i < numRows - 1; i++) {
0328: result += ldiagRe[i] * ldiagRe[i] + ldiagIm[i] * ldiagIm[i]
0329: + diagRe[i] * diagRe[i] + diagIm[i] * diagIm[i]
0330: + udiagRe[i] * udiagRe[i] + udiagIm[i] * udiagIm[i];
0331: }
0332: result += ldiagRe[i] * ldiagRe[i] + ldiagIm[i] * ldiagIm[i]
0333: + diagRe[i] * diagRe[i] + diagIm[i] * diagIm[i];
0334: return Math.sqrt(result);
0335: }
0336:
0337: /**
0338: * Returns the operator norm.
0339: * @exception MaximumIterationsExceededException If it takes more than 50 iterations to determine an eigenvalue.
0340: */
0341: public double operatorNorm()
0342: throws MaximumIterationsExceededException {
0343: return Math
0344: .sqrt(ArrayMath
0345: .max(LinearMath
0346: .eigenvalueSolveHermitian((ComplexTridiagonalMatrix) (this
0347: .hermitianAdjoint()
0348: .multiply(this )))));
0349: }
0350:
0351: //============
0352: // OPERATIONS
0353: //============
0354:
0355: // ADDITION
0356:
0357: /**
0358: * Returns the addition of this matrix and another.
0359: * @param m a complex matrix
0360: * @exception MatrixDimensionException If the matrices are different sizes.
0361: */
0362: public AbstractComplexSquareMatrix add(
0363: final AbstractComplexSquareMatrix m) {
0364: if (m instanceof ComplexTridiagonalMatrix)
0365: return add((ComplexTridiagonalMatrix) m);
0366: if (m instanceof TridiagonalMatrix)
0367: return addTridiagonal(m);
0368: if (m instanceof ComplexSquareMatrix)
0369: return add((ComplexSquareMatrix) m);
0370:
0371: if (numRows == m.rows() && numCols == m.columns()) {
0372: final double arrayRe[][] = new double[numRows][numCols];
0373: final double arrayIm[][] = new double[numRows][numCols];
0374: for (int i = 0; i < numRows; i++) {
0375: Complex elem = getElement(i, 0).add(m.getElement(i, 0));
0376: arrayRe[i][0] = elem.real();
0377: arrayIm[i][0] = elem.imag();
0378: for (int j = 1; j < numCols; j++) {
0379: elem = getElement(i, j).add(m.getElement(i, j));
0380: arrayRe[i][j] = elem.real();
0381: arrayIm[i][j] = elem.imag();
0382: }
0383: }
0384: return new ComplexSquareMatrix(arrayRe, arrayIm);
0385: } else {
0386: throw new MatrixDimensionException(
0387: "Matrices are different sizes.");
0388: }
0389: }
0390:
0391: public ComplexSquareMatrix add(final ComplexSquareMatrix m) {
0392: if (numRows == m.numRows && numCols == m.numCols) {
0393: final double arrayRe[][] = new double[numRows][numRows];
0394: final double arrayIm[][] = new double[numRows][numRows];
0395: for (int i = 0; i < numRows; i++) {
0396: System.arraycopy(m.matrixRe[i], 0, arrayRe[i], 0,
0397: numRows);
0398: System.arraycopy(m.matrixIm[i], 0, arrayIm[i], 0,
0399: numRows);
0400: }
0401: arrayRe[0][0] += diagRe[0];
0402: arrayIm[0][0] += diagIm[0];
0403: arrayRe[0][1] += udiagRe[0];
0404: arrayIm[0][1] += udiagIm[0];
0405: int n = numCols - 1;
0406: for (int i = 1; i < n; i++) {
0407: arrayRe[i][i - 1] += ldiagRe[i];
0408: arrayIm[i][i - 1] += ldiagIm[i];
0409: arrayRe[i][i] += diagRe[i];
0410: arrayIm[i][i] += diagIm[i];
0411: arrayRe[i][i + 1] += udiagRe[i];
0412: arrayIm[i][i + 1] += udiagIm[i];
0413: }
0414: arrayRe[n][n - 1] += ldiagRe[n];
0415: arrayIm[n][n - 1] += ldiagIm[n];
0416: arrayRe[n][n] += diagRe[n];
0417: arrayIm[n][n] += diagIm[n];
0418: return new ComplexSquareMatrix(arrayRe, arrayIm);
0419: } else
0420: throw new MatrixDimensionException(
0421: "Matrices are different sizes.");
0422: }
0423:
0424: public ComplexTridiagonalMatrix add(final ComplexTridiagonalMatrix m) {
0425: int mRow = numRows;
0426: if (mRow == m.numRows) {
0427: final ComplexTridiagonalMatrix ans = new ComplexTridiagonalMatrix(
0428: mRow);
0429: ans.diagRe[0] = diagRe[0] + m.diagRe[0];
0430: ans.diagIm[0] = diagIm[0] + m.diagIm[0];
0431: ans.udiagRe[0] = udiagRe[0] + m.udiagRe[0];
0432: ans.udiagIm[0] = udiagIm[0] + m.udiagIm[0];
0433: mRow--;
0434: for (int i = 1; i < mRow; i++) {
0435: ans.ldiagRe[i] = ldiagRe[i] + m.ldiagRe[i];
0436: ans.ldiagIm[i] = ldiagIm[i] + m.ldiagIm[i];
0437: ans.diagRe[i] = diagRe[i] + m.diagRe[i];
0438: ans.diagIm[i] = diagIm[i] + m.diagIm[i];
0439: ans.udiagRe[i] = udiagRe[i] + m.udiagRe[i];
0440: ans.udiagIm[i] = udiagIm[i] + m.udiagIm[i];
0441: }
0442: ans.ldiagRe[mRow] = ldiagRe[mRow] + m.ldiagRe[mRow];
0443: ans.ldiagIm[mRow] = ldiagIm[mRow] + m.ldiagIm[mRow];
0444: ans.diagRe[mRow] = diagRe[mRow] + m.diagRe[mRow];
0445: ans.diagIm[mRow] = diagIm[mRow] + m.diagIm[mRow];
0446: return ans;
0447: } else
0448: throw new MatrixDimensionException(
0449: "Matrices are different sizes.");
0450: }
0451:
0452: private ComplexTridiagonalMatrix addTridiagonal(
0453: final AbstractComplexSquareMatrix m) {
0454: int mRow = numRows;
0455: if (mRow == m.rows()) {
0456: final ComplexTridiagonalMatrix ans = new ComplexTridiagonalMatrix(
0457: mRow);
0458: Complex elem = m.getElement(0, 0);
0459: ans.diagRe[0] = diagRe[0] + elem.real();
0460: ans.diagIm[0] = diagIm[0] + elem.imag();
0461: elem = m.getElement(0, 1);
0462: ans.udiagRe[0] = udiagRe[0] + elem.real();
0463: ans.udiagIm[0] = udiagIm[0] + elem.imag();
0464: mRow--;
0465: for (int i = 1; i < mRow; i++) {
0466: elem = m.getElement(i, i - 1);
0467: ans.ldiagRe[i] = ldiagRe[i] + elem.real();
0468: ans.ldiagIm[i] = ldiagIm[i] + elem.imag();
0469: elem = m.getElement(i, i);
0470: ans.diagRe[i] = diagRe[i] + elem.real();
0471: ans.diagIm[i] = diagIm[i] + elem.imag();
0472: elem = m.getElement(i, i + 1);
0473: ans.udiagRe[i] = udiagRe[i] + elem.real();
0474: ans.udiagIm[i] = udiagIm[i] + elem.imag();
0475: }
0476: elem = m.getElement(mRow, mRow - 1);
0477: ans.ldiagRe[mRow] = ldiagRe[mRow] + elem.real();
0478: ans.ldiagIm[mRow] = ldiagIm[mRow] + elem.imag();
0479: elem = m.getElement(mRow, mRow);
0480: ans.diagRe[mRow] = diagRe[mRow] + elem.real();
0481: ans.diagIm[mRow] = diagIm[mRow] + elem.imag();
0482: return ans;
0483: } else {
0484: throw new MatrixDimensionException(
0485: "Matrices are different sizes.");
0486: }
0487: }
0488:
0489: // SUBTRACTION
0490:
0491: /**
0492: * Returns the subtraction of this matrix and another.
0493: * @param m a complex matrix
0494: * @exception MatrixDimensionException If the matrices are different sizes.
0495: */
0496: public AbstractComplexSquareMatrix subtract(
0497: final AbstractComplexSquareMatrix m) {
0498: if (m instanceof ComplexTridiagonalMatrix)
0499: return subtract((ComplexTridiagonalMatrix) m);
0500: if (m instanceof TridiagonalMatrix)
0501: return subtractTridiagonal(m);
0502: if (m instanceof ComplexSquareMatrix)
0503: return subtract((ComplexSquareMatrix) m);
0504:
0505: if (numRows == m.rows() && numCols == m.columns()) {
0506: final double arrayRe[][] = new double[numRows][numCols];
0507: final double arrayIm[][] = new double[numRows][numCols];
0508: Complex elem;
0509: for (int i = 0; i < numRows; i++) {
0510: elem = getElement(i, 0).subtract(m.getElement(i, 0));
0511: arrayRe[i][0] = elem.real();
0512: arrayIm[i][0] = elem.imag();
0513: for (int j = 1; j < numCols; j++) {
0514: elem = getElement(i, j)
0515: .subtract(m.getElement(i, j));
0516: arrayRe[i][j] = elem.real();
0517: arrayIm[i][j] = elem.imag();
0518: }
0519: }
0520: return new ComplexSquareMatrix(arrayRe, arrayIm);
0521: } else {
0522: throw new MatrixDimensionException(
0523: "Matrices are different sizes.");
0524: }
0525: }
0526:
0527: public ComplexSquareMatrix subtract(final ComplexSquareMatrix m) {
0528: if (numRows == m.numRows && numCols == m.numCols) {
0529: final double arrayRe[][] = new double[numRows][numCols];
0530: final double arrayIm[][] = new double[numRows][numCols];
0531: Complex elem;
0532: for (int j, i = 0; i < numRows; i++) {
0533: elem = getElement(i, 0);
0534: arrayRe[i][0] = elem.real() - m.matrixRe[i][0];
0535: arrayIm[i][0] = elem.imag() - m.matrixIm[i][0];
0536: for (j = 1; j < numCols; j++) {
0537: elem = getElement(i, j);
0538: arrayRe[i][j] = elem.real() - m.matrixRe[i][j];
0539: arrayIm[i][j] = elem.imag() - m.matrixIm[i][j];
0540: }
0541: }
0542: return new ComplexSquareMatrix(arrayRe, arrayIm);
0543: } else
0544: throw new MatrixDimensionException(
0545: "Matrices are different sizes.");
0546: }
0547:
0548: /**
0549: * Returns the subtraction of this matrix and another.
0550: * @param m a complex tridiagonal matrix
0551: * @exception MatrixDimensionException If the matrices are different sizes.
0552: */
0553: public ComplexTridiagonalMatrix subtract(
0554: final ComplexTridiagonalMatrix m) {
0555: int mRow = numRows;
0556: if (mRow == m.numRows) {
0557: final ComplexTridiagonalMatrix ans = new ComplexTridiagonalMatrix(
0558: mRow);
0559: ans.diagRe[0] = diagRe[0] - m.diagRe[0];
0560: ans.diagIm[0] = diagIm[0] - m.diagIm[0];
0561: ans.udiagRe[0] = udiagRe[0] - m.udiagRe[0];
0562: ans.udiagIm[0] = udiagIm[0] - m.udiagIm[0];
0563: mRow--;
0564: for (int i = 1; i < mRow; i++) {
0565: ans.ldiagRe[i] = ldiagRe[i] - m.ldiagRe[i];
0566: ans.ldiagIm[i] = ldiagIm[i] - m.ldiagIm[i];
0567: ans.diagRe[i] = diagRe[i] - m.diagRe[i];
0568: ans.diagIm[i] = diagIm[i] - m.diagIm[i];
0569: ans.udiagRe[i] = udiagRe[i] - m.udiagRe[i];
0570: ans.udiagIm[i] = udiagIm[i] - m.udiagIm[i];
0571: }
0572: ans.ldiagRe[mRow] = ldiagRe[mRow] - m.ldiagRe[mRow];
0573: ans.ldiagIm[mRow] = ldiagIm[mRow] - m.ldiagIm[mRow];
0574: ans.diagRe[mRow] = diagRe[mRow] - m.diagRe[mRow];
0575: ans.diagIm[mRow] = diagIm[mRow] - m.diagIm[mRow];
0576: return ans;
0577: } else
0578: throw new MatrixDimensionException(
0579: "Matrices are different sizes.");
0580: }
0581:
0582: private ComplexTridiagonalMatrix subtractTridiagonal(
0583: final AbstractComplexSquareMatrix m) {
0584: int mRow = numRows;
0585: if (mRow == m.rows()) {
0586: final ComplexTridiagonalMatrix ans = new ComplexTridiagonalMatrix(
0587: mRow);
0588: Complex elem = m.getElement(0, 0);
0589: ans.diagRe[0] = diagRe[0] - elem.real();
0590: ans.diagIm[0] = diagIm[0] - elem.imag();
0591: elem = m.getElement(0, 1);
0592: ans.udiagRe[0] = udiagRe[0] - elem.real();
0593: ans.udiagIm[0] = udiagIm[0] - elem.imag();
0594: mRow--;
0595: for (int i = 1; i < mRow; i++) {
0596: elem = m.getElement(i, i - 1);
0597: ans.ldiagRe[i] = ldiagRe[i] - elem.real();
0598: ans.ldiagIm[i] = ldiagIm[i] - elem.imag();
0599: elem = m.getElement(i, i);
0600: ans.diagRe[i] = diagRe[i] - elem.real();
0601: ans.diagIm[i] = diagIm[i] - elem.imag();
0602: elem = m.getElement(i, i + 1);
0603: ans.udiagRe[i] = udiagRe[i] - elem.real();
0604: ans.udiagIm[i] = udiagIm[i] - elem.imag();
0605: }
0606: elem = m.getElement(mRow, mRow - 1);
0607: ans.ldiagRe[mRow] = ldiagRe[mRow] - elem.real();
0608: ans.ldiagIm[mRow] = ldiagIm[mRow] - elem.imag();
0609: elem = m.getElement(mRow, mRow);
0610: ans.diagRe[mRow] = diagRe[mRow] - elem.real();
0611: ans.diagIm[mRow] = diagIm[mRow] - elem.imag();
0612: return ans;
0613: } else {
0614: throw new MatrixDimensionException(
0615: "Matrices are different sizes.");
0616: }
0617: }
0618:
0619: // SCALAR MULTIPLICATION
0620:
0621: /**
0622: * Returns the multiplication of this matrix by a scalar.
0623: * @param z a complex number
0624: * @return a complex tridiagonal matrix
0625: */
0626: public AbstractComplexMatrix scalarMultiply(final Complex z) {
0627: final double real = z.real();
0628: final double imag = z.imag();
0629: int mRow = numRows;
0630: final ComplexTridiagonalMatrix ans = new ComplexTridiagonalMatrix(
0631: mRow);
0632: ans.diagRe[0] = real * diagRe[0] - imag * diagIm[0];
0633: ans.diagIm[0] = imag * diagRe[0] + real * diagIm[0];
0634: ans.udiagRe[0] = real * udiagRe[0] - imag * udiagIm[0];
0635: ans.udiagIm[0] = imag * udiagRe[0] + real * udiagIm[0];
0636: mRow--;
0637: for (int i = 1; i < mRow; i++) {
0638: ans.ldiagRe[i] = real * ldiagRe[i] - imag * ldiagIm[i];
0639: ans.ldiagIm[i] = imag * ldiagRe[i] + real * ldiagIm[i];
0640: ans.diagRe[i] = real * diagRe[i] - imag * diagIm[i];
0641: ans.diagIm[i] = imag * diagRe[i] + real * diagIm[i];
0642: ans.udiagRe[i] = real * udiagRe[i] - imag * udiagIm[i];
0643: ans.udiagIm[i] = imag * udiagRe[i] + real * udiagIm[i];
0644: }
0645: ans.ldiagRe[mRow] = real * ldiagRe[mRow] - imag * ldiagIm[mRow];
0646: ans.ldiagIm[mRow] = imag * ldiagRe[mRow] + real * ldiagIm[mRow];
0647: ans.diagRe[mRow] = real * diagRe[mRow] - imag * diagIm[mRow];
0648: ans.diagIm[mRow] = imag * diagRe[mRow] + real * diagIm[mRow];
0649: return ans;
0650: }
0651:
0652: /**
0653: * Returns the multiplication of this matrix by a scalar.
0654: * @param x a double
0655: * @return a complex tridiagonal matrix
0656: */
0657: public AbstractComplexMatrix scalarMultiply(final double x) {
0658: int mRow = numRows;
0659: final ComplexTridiagonalMatrix ans = new ComplexTridiagonalMatrix(
0660: mRow);
0661: ans.diagRe[0] = x * diagRe[0];
0662: ans.diagIm[0] = x * diagIm[0];
0663: ans.udiagRe[0] = x * udiagRe[0];
0664: ans.udiagIm[0] = x * udiagIm[0];
0665: mRow--;
0666: for (int i = 1; i < mRow; i++) {
0667: ans.ldiagRe[i] = x * ldiagRe[i];
0668: ans.ldiagIm[i] = x * ldiagIm[i];
0669: ans.diagRe[i] = x * diagRe[i];
0670: ans.diagIm[i] = x * diagIm[i];
0671: ans.udiagRe[i] = x * udiagRe[i];
0672: ans.udiagIm[i] = x * udiagIm[i];
0673: }
0674: ans.ldiagRe[mRow] = x * ldiagRe[mRow];
0675: ans.ldiagIm[mRow] = x * ldiagIm[mRow];
0676: ans.diagRe[mRow] = x * diagRe[mRow];
0677: ans.diagIm[mRow] = x * diagIm[mRow];
0678: return ans;
0679: }
0680:
0681: // MATRIX MULTIPLICATION
0682:
0683: /**
0684: * Returns the multiplication of a vector by this matrix.
0685: * @param v a complex vector
0686: * @exception DimensionException If the matrix and vector are incompatible.
0687: */
0688: public AbstractComplexVector multiply(final AbstractComplexVector v) {
0689: int mRow = numRows;
0690: if (mRow == v.dimension()) {
0691: final double arrayRe[] = new double[mRow];
0692: final double arrayIm[] = new double[mRow];
0693: Complex comp;
0694: comp = v.getComponent(0);
0695: arrayRe[0] = (diagRe[0] * comp.real() - diagIm[0]
0696: * comp.imag());
0697: arrayIm[0] = (diagIm[0] * comp.real() + diagRe[0]
0698: * comp.imag());
0699: comp = v.getComponent(1);
0700: arrayRe[0] += (udiagRe[0] * comp.real() - udiagIm[0]
0701: * comp.imag());
0702: arrayIm[0] += (udiagIm[0] * comp.real() + udiagRe[0]
0703: * comp.imag());
0704: mRow--;
0705: for (int i = 1; i < mRow; i++) {
0706: comp = v.getComponent(i - 1);
0707: arrayRe[i] = (ldiagRe[i] * comp.real() - ldiagIm[i]
0708: * comp.imag());
0709: arrayIm[i] = (ldiagIm[i] * comp.real() + ldiagRe[i]
0710: * comp.imag());
0711: comp = v.getComponent(i);
0712: arrayRe[i] += (diagRe[i] * comp.real() - diagIm[i]
0713: * comp.imag());
0714: arrayIm[i] += (diagIm[i] * comp.real() + diagRe[i]
0715: * comp.imag());
0716: comp = v.getComponent(i + 1);
0717: arrayRe[i] += (udiagRe[i] * comp.real() - udiagIm[i]
0718: * comp.imag());
0719: arrayIm[i] += (udiagIm[i] * comp.real() + udiagRe[i]
0720: * comp.imag());
0721: }
0722: comp = v.getComponent(mRow - 1);
0723: arrayRe[mRow] = (ldiagRe[mRow] * comp.real() - ldiagIm[mRow]
0724: * comp.imag());
0725: arrayIm[mRow] = (ldiagIm[mRow] * comp.real() + ldiagRe[mRow]
0726: * comp.imag());
0727: comp = v.getComponent(mRow);
0728: arrayRe[mRow] += (diagRe[mRow] * comp.real() - diagIm[mRow]
0729: * comp.imag());
0730: arrayIm[mRow] += (diagIm[mRow] * comp.real() + diagRe[mRow]
0731: * comp.imag());
0732: return new ComplexVector(arrayRe, arrayIm);
0733: } else
0734: throw new DimensionException(
0735: "Matrix and vector are incompatible.");
0736: }
0737:
0738: /**
0739: * Returns the multiplication of this matrix and another.
0740: * @param m a complex matrix
0741: * @exception MatrixDimensionException If the matrices are different sizes.
0742: */
0743: public AbstractComplexSquareMatrix multiply(
0744: final AbstractComplexSquareMatrix m) {
0745: if (m instanceof ComplexTridiagonalMatrix)
0746: return multiply((ComplexTridiagonalMatrix) m);
0747: if (m instanceof TridiagonalMatrix)
0748: return multiplyTridiagonal(m);
0749: if (m instanceof ComplexSquareMatrix)
0750: return multiply((ComplexSquareMatrix) m);
0751:
0752: if (numCols == m.rows()) {
0753: final int mColumns = m.columns();
0754: final double arrayRe[][] = new double[numRows][mColumns];
0755: final double arrayIm[][] = new double[numRows][mColumns];
0756: final int lastRow = numRows - 1;
0757: for (int j = 0; j < numRows; j++) {
0758: arrayRe[0][j] = (diagRe[0] * m.getRealElement(0, j) - diagIm[0]
0759: * m.getImagElement(0, j))
0760: + (udiagRe[0] * m.getRealElement(1, j) - udiagIm[0]
0761: * m.getImagElement(1, j));
0762: arrayIm[0][j] = (diagIm[0] * m.getRealElement(0, j) + diagRe[0]
0763: * m.getImagElement(0, j))
0764: + (udiagIm[0] * m.getRealElement(1, j) + udiagRe[0]
0765: * m.getImagElement(1, j));
0766: for (int i = 1; i < lastRow; i++) {
0767: arrayRe[i][j] = (ldiagRe[i]
0768: * m.getRealElement(i - 1, j) - ldiagIm[i]
0769: * m.getImagElement(i - 1, j))
0770: + (diagRe[i] * m.getRealElement(i, j) - diagIm[i]
0771: * m.getImagElement(i, j))
0772: + (udiagRe[i] * m.getRealElement(i + 1, j) - udiagIm[i]
0773: * m.getImagElement(i + 1, j));
0774: arrayIm[i][j] = (ldiagIm[i]
0775: * m.getRealElement(i - 1, j) + ldiagRe[i]
0776: * m.getImagElement(i - 1, j))
0777: + (diagIm[i] * m.getRealElement(i, j) + diagRe[i]
0778: * m.getImagElement(i, j))
0779: + (udiagIm[i] * m.getRealElement(i + 1, j) + udiagRe[i]
0780: * m.getImagElement(i + 1, j));
0781: }
0782: arrayRe[lastRow][j] = (ldiagRe[lastRow]
0783: * m.getRealElement(lastRow - 1, j) - ldiagIm[lastRow]
0784: * m.getImagElement(lastRow - 1, j))
0785: + (diagRe[lastRow]
0786: * m.getRealElement(lastRow, j) - diagIm[lastRow]
0787: * m.getImagElement(lastRow, j));
0788: arrayIm[lastRow][j] = (ldiagIm[lastRow]
0789: * m.getRealElement(lastRow - 1, j) + ldiagRe[lastRow]
0790: * m.getImagElement(lastRow - 1, j))
0791: + (diagIm[lastRow]
0792: * m.getRealElement(lastRow, j) + diagRe[lastRow]
0793: * m.getImagElement(lastRow, j));
0794: }
0795: return new ComplexSquareMatrix(arrayRe, arrayIm);
0796: } else {
0797: throw new MatrixDimensionException("Incompatible matrices.");
0798: }
0799: }
0800:
0801: public ComplexSquareMatrix multiply(final ComplexSquareMatrix m) {
0802: if (numCols == m.numRows) {
0803: final double arrayRe[][] = new double[numRows][m.numCols];
0804: final double arrayIm[][] = new double[numRows][m.numCols];
0805: final int lastRow = numRows - 1;
0806: for (int j = 0; j < numRows; j++) {
0807: arrayRe[0][j] = (diagRe[0] * m.matrixRe[0][j] - diagIm[0]
0808: * m.matrixIm[0][j])
0809: + (udiagRe[0] * m.matrixRe[1][j] - udiagIm[0]
0810: * m.matrixIm[1][j]);
0811: arrayIm[0][j] = (diagIm[0] * m.matrixRe[0][j] + diagRe[0]
0812: * m.matrixIm[0][j])
0813: + (udiagIm[0] * m.matrixRe[1][j] + udiagRe[0]
0814: * m.matrixIm[1][j]);
0815: for (int i = 1; i < lastRow; i++) {
0816: arrayRe[i][j] = (ldiagRe[i] * m.matrixRe[i - 1][j] - ldiagIm[i]
0817: * m.matrixIm[i - 1][j])
0818: + (diagRe[i] * m.matrixRe[i][j] - diagIm[i]
0819: * m.matrixIm[i][j])
0820: + (udiagRe[i] * m.matrixRe[i + 1][j] - udiagIm[i]
0821: * m.matrixIm[i + 1][j]);
0822: arrayIm[i][j] = (ldiagIm[i] * m.matrixRe[i - 1][j] + ldiagRe[i]
0823: * m.matrixIm[i - 1][j])
0824: + (diagIm[i] * m.matrixRe[i][j] + diagRe[i]
0825: * m.matrixIm[i][j])
0826: + (udiagIm[i] * m.matrixRe[i + 1][j] + udiagRe[i]
0827: * m.matrixIm[i + 1][j]);
0828: }
0829: arrayRe[lastRow][j] = (ldiagRe[lastRow]
0830: * m.matrixRe[lastRow - 1][j] - ldiagIm[lastRow]
0831: * m.matrixIm[lastRow - 1][j])
0832: + (diagRe[lastRow] * m.matrixRe[lastRow][j] - diagIm[lastRow]
0833: * m.matrixIm[lastRow][j]);
0834: arrayIm[lastRow][j] = (ldiagIm[lastRow]
0835: * m.matrixRe[lastRow - 1][j] + ldiagRe[lastRow]
0836: * m.matrixIm[lastRow - 1][j])
0837: + (diagIm[lastRow] * m.matrixRe[lastRow][j] + diagRe[lastRow]
0838: * m.matrixIm[lastRow][j]);
0839: }
0840: return new ComplexSquareMatrix(arrayRe, arrayIm);
0841: } else {
0842: throw new MatrixDimensionException("Incompatible matrices.");
0843: }
0844: }
0845:
0846: /**
0847: * Returns the multiplication of this matrix and another.
0848: * @param m a complex tridiagonal matrix
0849: * @exception MatrixDimensionException If the matrices are different sizes.
0850: */
0851: public ComplexSquareMatrix multiply(final ComplexTridiagonalMatrix m) {
0852: int mRow = numRows;
0853: if (numCols == m.numRows) {
0854: final double arrayRe[][] = new double[mRow][mRow];
0855: final double arrayIm[][] = new double[mRow][mRow];
0856: arrayRe[0][0] = (diagRe[0] * m.diagRe[0] - diagIm[0]
0857: * m.diagIm[0])
0858: + (udiagRe[0] * m.ldiagRe[1] - udiagIm[0]
0859: * m.ldiagIm[1]);
0860: arrayIm[0][0] = (diagIm[0] * m.diagRe[0] + diagRe[0]
0861: * m.diagIm[0])
0862: + (udiagIm[0] * m.ldiagRe[1] + udiagRe[0]
0863: * m.ldiagIm[1]);
0864: arrayRe[0][1] = (diagRe[0] * m.udiagRe[0] - diagIm[0]
0865: * m.udiagIm[0])
0866: + (udiagRe[0] * m.diagRe[1] - udiagIm[0]
0867: * m.diagIm[1]);
0868: arrayIm[0][1] = (diagIm[0] * m.udiagRe[0] + diagRe[0]
0869: * m.udiagIm[0])
0870: + (udiagIm[0] * m.diagRe[1] + udiagRe[0]
0871: * m.diagIm[1]);
0872: arrayRe[0][2] = (udiagRe[0] * m.udiagRe[1] - udiagIm[0]
0873: * m.udiagIm[1]);
0874: arrayIm[0][2] = (udiagIm[0] * m.udiagRe[1] + udiagRe[0]
0875: * m.udiagIm[1]);
0876: if (mRow > 3) {
0877: arrayRe[1][0] = (ldiagRe[1] * m.diagRe[0] - ldiagIm[1]
0878: * m.diagIm[0])
0879: + (diagRe[1] * m.ldiagRe[1] - diagIm[1]
0880: * m.ldiagIm[1]);
0881: arrayIm[1][0] = (ldiagIm[1] * m.diagRe[0] + ldiagRe[1]
0882: * m.diagIm[0])
0883: + (diagIm[1] * m.ldiagRe[1] + diagRe[1]
0884: * m.ldiagIm[1]);
0885: arrayRe[1][1] = (ldiagRe[1] * m.udiagRe[0] - ldiagIm[1]
0886: * m.udiagIm[0])
0887: + (diagRe[1] * m.diagRe[1] - diagIm[1]
0888: * m.diagIm[1])
0889: + (udiagRe[1] * m.ldiagRe[2] - udiagIm[1]
0890: * m.ldiagIm[2]);
0891: arrayIm[1][1] = (ldiagIm[1] * m.udiagRe[0] + ldiagRe[1]
0892: * m.udiagIm[0])
0893: + (diagIm[1] * m.diagRe[1] + diagRe[1]
0894: * m.diagIm[1])
0895: + (udiagIm[1] * m.ldiagRe[2] + udiagRe[1]
0896: * m.ldiagIm[2]);
0897: arrayRe[1][2] = (diagRe[1] * m.udiagRe[1] - diagIm[1]
0898: * m.udiagIm[1])
0899: + (udiagRe[1] * m.diagRe[2] - udiagIm[1]
0900: * m.diagIm[2]);
0901: arrayIm[1][2] = (diagIm[1] * m.udiagRe[1] + diagRe[1]
0902: * m.udiagIm[1])
0903: + (udiagIm[1] * m.diagRe[2] + udiagRe[1]
0904: * m.diagIm[2]);
0905: arrayRe[1][3] = (udiagRe[1] * m.udiagRe[2] - udiagIm[1]
0906: * m.udiagIm[2]);
0907: arrayIm[1][3] = (udiagIm[1] * m.udiagRe[2] + udiagRe[1]
0908: * m.udiagIm[2]);
0909: }
0910: if (mRow == 3) {
0911: arrayRe[1][0] = (ldiagRe[1] * m.diagRe[0] - ldiagIm[1]
0912: * m.diagIm[0])
0913: + (diagRe[1] * m.ldiagRe[1] - diagIm[1]
0914: * m.ldiagIm[1]);
0915: arrayIm[1][0] = (ldiagIm[1] * m.diagRe[0] + ldiagRe[1]
0916: * m.diagIm[0])
0917: + (diagIm[1] * m.ldiagRe[1] + diagRe[1]
0918: * m.ldiagIm[1]);
0919: arrayRe[1][1] = (ldiagRe[1] * m.udiagRe[0] - ldiagIm[1]
0920: * m.udiagIm[0])
0921: + (diagRe[1] * m.diagRe[1] - diagIm[1]
0922: * m.diagIm[1])
0923: + (udiagRe[1] * m.ldiagRe[2] - udiagIm[1]
0924: * m.ldiagIm[2]);
0925: arrayIm[1][1] = (ldiagIm[1] * m.udiagRe[0] + ldiagRe[1]
0926: * m.udiagIm[0])
0927: + (diagIm[1] * m.diagRe[1] + diagRe[1]
0928: * m.diagIm[1])
0929: + (udiagIm[1] * m.ldiagRe[2] + udiagRe[1]
0930: * m.ldiagIm[2]);
0931: arrayRe[1][2] = (diagRe[1] * m.udiagRe[1] - diagIm[1]
0932: * m.udiagIm[1])
0933: + (udiagRe[1] * m.diagRe[2] - udiagIm[1]
0934: * m.diagIm[2]);
0935: arrayIm[1][2] = (diagIm[1] * m.udiagRe[1] + diagRe[1]
0936: * m.udiagIm[1])
0937: + (udiagIm[1] * m.diagRe[2] + udiagRe[1]
0938: * m.diagIm[2]);
0939: } else if (mRow > 4) {
0940: for (int i = 2; i < mRow - 2; i++) {
0941: arrayRe[i][i - 2] = (ldiagRe[i] * m.ldiagRe[i - 1] - ldiagIm[i]
0942: * m.ldiagIm[i - 1]);
0943: arrayIm[i][i - 2] = (ldiagIm[i] * m.ldiagRe[i - 1] + ldiagRe[i]
0944: * m.ldiagIm[i - 1]);
0945: arrayRe[i][i - 1] = (ldiagRe[i] * m.diagRe[i - 1] - ldiagIm[i]
0946: * m.diagIm[i - 1])
0947: + (diagRe[i] * m.ldiagRe[i] - diagIm[i]
0948: * m.ldiagIm[i]);
0949: arrayIm[i][i - 1] = (ldiagIm[i] * m.diagRe[i - 1] + ldiagRe[i]
0950: * m.diagIm[i - 1])
0951: + (diagIm[i] * m.ldiagRe[i] + diagRe[i]
0952: * m.ldiagIm[i]);
0953: arrayRe[i][i] = (ldiagRe[i] * m.udiagRe[i - 1] - ldiagIm[i]
0954: * m.udiagIm[i - 1])
0955: + (diagRe[i] * m.diagRe[i] - diagIm[i]
0956: * m.diagIm[i])
0957: + (udiagRe[i] * m.ldiagRe[i + 1] - udiagIm[i]
0958: * m.ldiagIm[i + 1]);
0959: arrayIm[i][i] = (ldiagIm[i] * m.udiagRe[i - 1] + ldiagRe[i]
0960: * m.udiagIm[i - 1])
0961: + (diagIm[i] * m.diagRe[i] + diagRe[i]
0962: * m.diagIm[i])
0963: + (udiagIm[i] * m.ldiagRe[i + 1] + udiagRe[i]
0964: * m.ldiagIm[i + 1]);
0965: arrayRe[i][i + 1] = (diagRe[i] * m.udiagRe[i] - diagIm[i]
0966: * m.udiagIm[i])
0967: + (udiagRe[i] * m.diagRe[i + 1] - udiagIm[i]
0968: * m.diagIm[i + 1]);
0969: arrayIm[i][i + 1] = (diagIm[i] * m.udiagRe[i] + diagRe[i]
0970: * m.udiagIm[i])
0971: + (udiagIm[i] * m.diagRe[i + 1] + udiagRe[i]
0972: * m.diagIm[i + 1]);
0973: arrayRe[i][i + 2] = (udiagRe[i] * m.udiagRe[i + 1] - udiagIm[i]
0974: * m.udiagIm[i + 1]);
0975: arrayIm[i][i + 2] = (udiagIm[i] * m.udiagRe[i + 1] + udiagRe[i]
0976: * m.udiagIm[i + 1]);
0977: }
0978: }
0979: if (mRow > 3) {
0980: arrayRe[mRow - 2][mRow - 4] = (ldiagRe[mRow - 2]
0981: * m.ldiagRe[mRow - 3] - ldiagIm[mRow - 2]
0982: * m.ldiagIm[mRow - 3]);
0983: arrayIm[mRow - 2][mRow - 4] = (ldiagIm[mRow - 2]
0984: * m.ldiagRe[mRow - 3] + ldiagRe[mRow - 2]
0985: * m.ldiagIm[mRow - 3]);
0986: arrayRe[mRow - 2][mRow - 3] = (ldiagRe[mRow - 2]
0987: * m.diagRe[mRow - 3] - ldiagIm[mRow - 2]
0988: * m.diagIm[mRow - 3])
0989: + (diagRe[mRow - 2] * m.ldiagRe[mRow - 2] - diagIm[mRow - 2]
0990: * m.ldiagIm[mRow - 2]);
0991: arrayIm[mRow - 2][mRow - 3] = (ldiagIm[mRow - 2]
0992: * m.diagRe[mRow - 3] + ldiagRe[mRow - 2]
0993: * m.diagIm[mRow - 3])
0994: + (diagIm[mRow - 2] * m.ldiagRe[mRow - 2] + diagRe[mRow - 2]
0995: * m.ldiagIm[mRow - 2]);
0996: arrayRe[mRow - 2][mRow - 2] = (ldiagRe[mRow - 2]
0997: * m.udiagRe[mRow - 3] - ldiagIm[mRow - 2]
0998: * m.udiagIm[mRow - 3])
0999: + (diagRe[mRow - 2] * m.diagRe[mRow - 2] - diagIm[mRow - 2]
1000: * m.diagIm[mRow - 2])
1001: + (udiagRe[mRow - 2] * m.ldiagRe[mRow - 1] - udiagIm[mRow - 2]
1002: * m.ldiagIm[mRow - 1]);
1003: arrayIm[mRow - 2][mRow - 2] = (ldiagIm[mRow - 2]
1004: * m.udiagRe[mRow - 3] + ldiagRe[mRow - 2]
1005: * m.udiagIm[mRow - 3])
1006: + (diagIm[mRow - 2] * m.diagRe[mRow - 2] + diagRe[mRow - 2]
1007: * m.diagIm[mRow - 2])
1008: + (udiagIm[mRow - 2] * m.ldiagRe[mRow - 1] + udiagRe[mRow - 2]
1009: * m.ldiagIm[mRow - 1]);
1010: arrayRe[mRow - 2][mRow - 1] = (diagRe[mRow - 2]
1011: * m.udiagRe[mRow - 2] - diagIm[mRow - 2]
1012: * m.udiagIm[mRow - 2])
1013: + (udiagRe[mRow - 2] * m.diagRe[mRow - 1] - udiagIm[mRow - 2]
1014: * m.diagIm[mRow - 1]);
1015: arrayIm[mRow - 2][mRow - 1] = (diagIm[mRow - 2]
1016: * m.udiagRe[mRow - 2] + diagRe[mRow - 2]
1017: * m.udiagIm[mRow - 2])
1018: + (udiagIm[mRow - 2] * m.diagRe[mRow - 1] + udiagRe[mRow - 2]
1019: * m.diagIm[mRow - 1]);
1020: }
1021: mRow--;
1022: arrayRe[mRow][mRow - 2] = (ldiagRe[mRow]
1023: * m.ldiagRe[mRow - 1] - ldiagIm[mRow]
1024: * m.ldiagIm[mRow - 1]);
1025: arrayIm[mRow][mRow - 2] = (ldiagIm[mRow]
1026: * m.ldiagRe[mRow - 1] + ldiagRe[mRow]
1027: * m.ldiagIm[mRow - 1]);
1028: arrayRe[mRow][mRow - 1] = (ldiagRe[mRow]
1029: * m.diagRe[mRow - 1] - ldiagIm[mRow]
1030: * m.diagIm[mRow - 1])
1031: + (diagRe[mRow] * m.ldiagRe[mRow] - diagIm[mRow]
1032: * m.ldiagIm[mRow]);
1033: arrayIm[mRow][mRow - 1] = (ldiagIm[mRow]
1034: * m.diagRe[mRow - 1] + ldiagRe[mRow]
1035: * m.diagIm[mRow - 1])
1036: + (diagIm[mRow] * m.ldiagRe[mRow] + diagRe[mRow]
1037: * m.ldiagIm[mRow]);
1038: arrayRe[mRow][mRow] = (ldiagRe[mRow] * m.udiagRe[mRow - 1] - ldiagIm[mRow]
1039: * m.udiagIm[mRow - 1])
1040: + (diagRe[mRow] * m.diagRe[mRow] - diagIm[mRow]
1041: * m.diagIm[mRow]);
1042: arrayIm[mRow][mRow] = (ldiagIm[mRow] * m.udiagRe[mRow - 1] + ldiagRe[mRow]
1043: * m.udiagIm[mRow - 1])
1044: + (diagIm[mRow] * m.diagRe[mRow] + diagRe[mRow]
1045: * m.diagIm[mRow]);
1046: return new ComplexSquareMatrix(arrayRe, arrayIm);
1047: } else
1048: throw new MatrixDimensionException("Incompatible matrices.");
1049: }
1050:
1051: private ComplexSquareMatrix multiplyTridiagonal(
1052: final AbstractComplexSquareMatrix m) {
1053: int mRow = numRows;
1054: if (numCols == m.rows()) {
1055: final double arrayRe[][] = new double[mRow][mRow];
1056: final double arrayIm[][] = new double[mRow][mRow];
1057: Complex elem1, elem2, elem3;
1058: elem1 = m.getElement(0, 0);
1059: elem2 = m.getElement(1, 0);
1060: arrayRe[0][0] = (diagRe[0] * elem1.real() - diagIm[0]
1061: * elem1.imag())
1062: + (udiagRe[0] * elem2.real() - udiagIm[0]
1063: * elem2.imag());
1064: arrayIm[0][0] = (diagIm[0] * elem1.real() + diagRe[0]
1065: * elem1.imag())
1066: + (udiagIm[0] * elem2.real() + udiagRe[0]
1067: * elem2.imag());
1068: elem1 = m.getElement(0, 1);
1069: elem2 = m.getElement(1, 1);
1070: arrayRe[0][1] = (diagRe[0] * elem1.real() - diagIm[0]
1071: * elem1.imag())
1072: + (udiagRe[0] * elem2.real() - udiagIm[0]
1073: * elem2.imag());
1074: arrayIm[0][1] = (diagIm[0] * elem1.real() + diagRe[0]
1075: * elem1.imag())
1076: + (udiagIm[0] * elem2.real() + udiagRe[0]
1077: * elem2.imag());
1078: elem1 = m.getElement(1, 2);
1079: arrayRe[0][2] = (udiagRe[0] * elem1.real() - udiagIm[0]
1080: * elem1.imag());
1081: arrayIm[0][2] = (udiagIm[0] * elem1.real() + udiagRe[0]
1082: * elem1.imag());
1083: if (mRow > 3) {
1084: elem1 = m.getElement(0, 0);
1085: elem2 = m.getElement(1, 0);
1086: arrayRe[1][0] = (ldiagRe[1] * elem1.real() - ldiagIm[1]
1087: * elem1.imag())
1088: + (diagRe[1] * elem2.real() - diagIm[1]
1089: * elem2.imag());
1090: arrayIm[1][0] = (ldiagIm[1] * elem1.real() + ldiagRe[1]
1091: * elem1.imag())
1092: + (diagIm[1] * elem2.real() + diagRe[1]
1093: * elem2.imag());
1094: elem1 = m.getElement(0, 1);
1095: elem2 = m.getElement(1, 1);
1096: elem3 = m.getElement(2, 1);
1097: arrayRe[1][1] = (ldiagRe[1] * elem1.real() - ldiagIm[1]
1098: * elem1.imag())
1099: + (diagRe[1] * elem2.real() - diagIm[1]
1100: * elem2.imag())
1101: + (udiagRe[1] * elem3.real() - udiagIm[1]
1102: * elem3.imag());
1103: arrayIm[1][1] = (ldiagIm[1] * elem1.real() + ldiagRe[1]
1104: * elem1.imag())
1105: + (diagIm[1] * elem2.real() + diagRe[1]
1106: * elem2.imag())
1107: + (udiagIm[1] * elem3.real() + udiagRe[1]
1108: * elem3.imag());
1109: elem1 = m.getElement(1, 2);
1110: elem2 = m.getElement(2, 2);
1111: arrayRe[1][2] = (diagRe[1] * elem1.real() - diagIm[1]
1112: * elem1.imag())
1113: + (udiagRe[1] * elem2.real() - udiagIm[1]
1114: * elem2.imag());
1115: arrayIm[1][2] = (diagIm[1] * elem1.real() + diagRe[1]
1116: * elem1.imag())
1117: + (udiagIm[1] * elem2.real() + udiagRe[1]
1118: * elem2.imag());
1119: elem1 = m.getElement(2, 3);
1120: arrayRe[1][3] = (udiagRe[1] * elem1.real() - udiagIm[1]
1121: * elem1.imag());
1122: arrayIm[1][3] = (udiagIm[1] * elem1.real() + udiagRe[1]
1123: * elem1.imag());
1124: }
1125: if (mRow == 3) {
1126: elem1 = m.getElement(0, 0);
1127: elem2 = m.getElement(1, 0);
1128: arrayRe[1][0] = (ldiagRe[1] * elem1.real() - ldiagIm[1]
1129: * elem1.imag())
1130: + (diagRe[1] * elem2.real() - diagIm[1]
1131: * elem2.imag());
1132: arrayIm[1][0] = (ldiagIm[1] * elem1.real() + ldiagRe[1]
1133: * elem1.imag())
1134: + (diagIm[1] * elem2.real() + diagRe[1]
1135: * elem2.imag());
1136: elem1 = m.getElement(0, 1);
1137: elem2 = m.getElement(1, 1);
1138: elem3 = m.getElement(2, 1);
1139: arrayRe[1][1] = (ldiagRe[1] * elem1.real() - ldiagIm[1]
1140: * elem1.imag())
1141: + (diagRe[1] * elem2.real() - diagIm[1]
1142: * elem2.imag())
1143: + (udiagRe[1] * elem3.real() - udiagIm[1]
1144: * elem3.imag());
1145: arrayIm[1][1] = (ldiagIm[1] * elem1.real() + ldiagRe[1]
1146: * elem1.imag())
1147: + (diagIm[1] * elem2.real() + diagRe[1]
1148: * elem2.imag())
1149: + (udiagIm[1] * elem3.real() + udiagRe[1]
1150: * elem3.imag());
1151: elem1 = m.getElement(1, 2);
1152: elem2 = m.getElement(2, 2);
1153: arrayRe[1][2] = (diagRe[1] * elem1.real() - diagIm[1]
1154: * elem1.imag())
1155: + (udiagRe[1] * elem2.real() - udiagIm[1]
1156: * elem2.imag());
1157: arrayIm[1][2] = (diagIm[1] * elem1.real() + diagRe[1]
1158: * elem1.imag())
1159: + (udiagIm[1] * elem2.real() + udiagRe[1]
1160: * elem2.imag());
1161: } else if (mRow > 4) {
1162: for (int i = 2; i < mRow - 2; i++) {
1163: elem1 = m.getElement(i - 1, i - 2);
1164: arrayRe[i][i - 2] = (ldiagRe[i] * elem1.real() - ldiagIm[i]
1165: * elem1.imag());
1166: arrayIm[i][i - 2] = (ldiagIm[i] * elem1.real() + ldiagRe[i]
1167: * elem1.imag());
1168: elem1 = m.getElement(i - 1, i - 1);
1169: elem2 = m.getElement(i, i - 1);
1170: arrayRe[i][i - 1] = (ldiagRe[i] * elem1.real() - ldiagIm[i]
1171: * elem1.imag())
1172: + (diagRe[i] * elem2.real() - diagIm[i]
1173: * elem2.imag());
1174: arrayIm[i][i - 1] = (ldiagIm[i] * elem1.real() + ldiagRe[i]
1175: * elem1.imag())
1176: + (diagIm[i] * elem2.real() + diagRe[i]
1177: * elem2.imag());
1178: elem1 = m.getElement(i - 1, i);
1179: elem2 = m.getElement(i, i);
1180: elem3 = m.getElement(i + 1, i);
1181: arrayRe[i][i] = (ldiagRe[i] * elem1.real() - ldiagIm[i]
1182: * elem1.imag())
1183: + (diagRe[i] * elem2.real() - diagIm[i]
1184: * elem2.imag())
1185: + (udiagRe[i] * elem3.real() - udiagIm[i]
1186: * elem3.imag());
1187: arrayIm[i][i] = (ldiagIm[i] * elem1.real() + ldiagRe[i]
1188: * elem1.imag())
1189: + (diagIm[i] * elem2.real() + diagRe[i]
1190: * elem2.imag())
1191: + (udiagIm[i] * elem3.real() + udiagRe[i]
1192: * elem3.imag());
1193: elem1 = m.getElement(i, i + 1);
1194: elem2 = m.getElement(i + 1, i + 1);
1195: arrayRe[i][i + 1] = (diagRe[i] * elem1.real() - diagIm[i]
1196: * elem1.imag())
1197: + (udiagRe[i] * elem2.real() - udiagIm[i]
1198: * elem2.imag());
1199: arrayIm[i][i + 1] = (diagIm[i] * elem1.real() + diagRe[i]
1200: * elem1.imag())
1201: + (udiagIm[i] * elem2.real() + udiagRe[i]
1202: * elem2.imag());
1203: elem1 = m.getElement(i + 1, i + 2);
1204: arrayRe[i][i + 2] = (udiagRe[i] * elem1.real() - udiagIm[i]
1205: * elem1.imag());
1206: arrayIm[i][i + 2] = (udiagIm[i] * elem1.real() + udiagRe[i]
1207: * elem1.imag());
1208: }
1209: }
1210: if (mRow > 3) {
1211: elem1 = m.getElement(mRow - 3, mRow - 4);
1212: arrayRe[mRow - 2][mRow - 4] = (ldiagRe[mRow - 2]
1213: * elem1.real() - ldiagIm[mRow - 2]
1214: * elem1.imag());
1215: arrayIm[mRow - 2][mRow - 4] = (ldiagIm[mRow - 2]
1216: * elem1.real() + ldiagRe[mRow - 2]
1217: * elem1.imag());
1218: elem1 = m.getElement(mRow - 3, mRow - 3);
1219: elem2 = m.getElement(mRow - 2, mRow - 3);
1220: arrayRe[mRow - 2][mRow - 3] = (ldiagRe[mRow - 2]
1221: * elem1.real() - ldiagIm[mRow - 2]
1222: * elem1.imag())
1223: + (diagRe[mRow - 2] * elem2.real() - diagIm[mRow - 2]
1224: * elem2.imag());
1225: arrayIm[mRow - 2][mRow - 3] = (ldiagIm[mRow - 2]
1226: * elem1.real() + ldiagRe[mRow - 2]
1227: * elem1.imag())
1228: + (diagIm[mRow - 2] * elem2.real() + diagRe[mRow - 2]
1229: * elem2.imag());
1230: elem1 = m.getElement(mRow - 3, mRow - 2);
1231: elem2 = m.getElement(mRow - 2, mRow - 2);
1232: elem3 = m.getElement(mRow - 1, mRow - 2);
1233: arrayRe[mRow - 2][mRow - 2] = (ldiagRe[mRow - 2]
1234: * elem1.real() - ldiagIm[mRow - 2]
1235: * elem1.imag())
1236: + (diagRe[mRow - 2] * elem2.real() - diagIm[mRow - 2]
1237: * elem2.imag())
1238: + (udiagRe[mRow - 2] * elem3.real() - udiagIm[mRow - 2]
1239: * elem3.imag());
1240: arrayIm[mRow - 2][mRow - 2] = (ldiagIm[mRow - 2]
1241: * elem1.real() + ldiagRe[mRow - 2]
1242: * elem1.imag())
1243: + (diagIm[mRow - 2] * elem2.real() + diagRe[mRow - 2]
1244: * elem2.imag())
1245: + (udiagIm[mRow - 2] * elem3.real() + udiagRe[mRow - 2]
1246: * elem3.imag());
1247: elem1 = m.getElement(mRow - 2, mRow - 1);
1248: elem2 = m.getElement(mRow - 1, mRow - 1);
1249: arrayRe[mRow - 2][mRow - 1] = (diagRe[mRow - 2]
1250: * elem1.real() - diagIm[mRow - 2]
1251: * elem1.imag())
1252: + (udiagRe[mRow - 2] * elem2.real() - udiagIm[mRow - 2]
1253: * elem2.imag());
1254: arrayIm[mRow - 2][mRow - 1] = (diagIm[mRow - 2]
1255: * elem1.real() + diagRe[mRow - 2]
1256: * elem1.imag())
1257: + (udiagIm[mRow - 2] * elem2.real() + udiagRe[mRow - 2]
1258: * elem2.imag());
1259: }
1260: mRow--;
1261: elem1 = m.getElement(mRow - 1, mRow - 2);
1262: arrayRe[mRow][mRow - 2] = (ldiagRe[mRow] * elem1.real() - ldiagIm[mRow]
1263: * elem1.imag());
1264: arrayIm[mRow][mRow - 2] = (ldiagIm[mRow] * elem1.real() + ldiagRe[mRow]
1265: * elem1.imag());
1266: elem1 = m.getElement(mRow - 1, mRow - 1);
1267: elem2 = m.getElement(mRow, mRow - 1);
1268: arrayRe[mRow][mRow - 1] = (ldiagRe[mRow] * elem1.real() - ldiagIm[mRow]
1269: * elem1.imag())
1270: + (diagRe[mRow] * elem2.real() - diagIm[mRow]
1271: * elem2.imag());
1272: arrayIm[mRow][mRow - 1] = (ldiagIm[mRow] * elem1.real() + ldiagRe[mRow]
1273: * elem1.imag())
1274: + (diagIm[mRow] * elem2.real() + diagRe[mRow]
1275: * elem2.imag());
1276: elem1 = m.getElement(mRow - 1, mRow);
1277: elem2 = m.getElement(mRow, mRow);
1278: arrayRe[mRow][mRow] = (ldiagRe[mRow] * elem1.real() - ldiagIm[mRow]
1279: * elem1.imag())
1280: + (diagRe[mRow] * elem2.real() - diagIm[mRow]
1281: * elem2.imag());
1282: arrayIm[mRow][mRow] = (ldiagIm[mRow] * elem1.real() + ldiagRe[mRow]
1283: * elem1.imag())
1284: + (diagIm[mRow] * elem2.real() + diagRe[mRow]
1285: * elem2.imag());
1286: return new ComplexSquareMatrix(arrayRe, arrayIm);
1287: } else {
1288: throw new MatrixDimensionException("Incompatible matrices.");
1289: }
1290: }
1291:
1292: // HERMITIAN ADJOINT
1293:
1294: /**
1295: * Returns the hermitian adjoint of this matrix.
1296: * @return a complex tridiagonal matrix
1297: */
1298: public AbstractComplexMatrix hermitianAdjoint() {
1299: int mRow = numRows;
1300: final ComplexTridiagonalMatrix ans = new ComplexTridiagonalMatrix(
1301: mRow);
1302: System
1303: .arraycopy(ldiagRe, 1, ans.udiagRe, 0,
1304: ldiagRe.length - 1);
1305: System.arraycopy(diagRe, 0, ans.diagRe, 0, diagRe.length);
1306: System
1307: .arraycopy(udiagRe, 0, ans.ldiagRe, 1,
1308: udiagRe.length - 1);
1309: ans.diagIm[0] = -diagIm[0];
1310: ans.ldiagIm[1] = -udiagIm[0];
1311: mRow--;
1312: for (int i = 1; i < mRow; i++) {
1313: ans.udiagIm[i - 1] = -ldiagIm[i];
1314: ans.diagIm[i] = -diagIm[i];
1315: ans.ldiagIm[i + 1] = -udiagIm[i];
1316: }
1317: ans.udiagIm[mRow - 1] = -ldiagIm[mRow];
1318: ans.diagIm[mRow] = -diagIm[mRow];
1319: return ans;
1320: }
1321:
1322: // CONJUGATE
1323:
1324: /**
1325: * Returns the complex conjugate of this matrix.
1326: * @return a complex tridiagonal matrix
1327: */
1328: public AbstractComplexMatrix conjugate() {
1329: int mRow = numRows;
1330: final ComplexTridiagonalMatrix ans = new ComplexTridiagonalMatrix(
1331: mRow);
1332: System
1333: .arraycopy(ldiagRe, 1, ans.ldiagRe, 0,
1334: ldiagRe.length - 1);
1335: System.arraycopy(diagRe, 0, ans.diagRe, 0, diagRe.length);
1336: System
1337: .arraycopy(udiagRe, 0, ans.udiagRe, 1,
1338: udiagRe.length - 1);
1339: ans.diagIm[0] = -diagIm[0];
1340: ans.udiagIm[0] = -udiagIm[0];
1341: mRow--;
1342: for (int i = 1; i < mRow; i++) {
1343: ans.ldiagIm[i] = -ldiagIm[i];
1344: ans.diagIm[i] = -diagIm[i];
1345: ans.udiagIm[i] = -udiagIm[i];
1346: }
1347: ans.ldiagIm[mRow] = -ldiagIm[mRow];
1348: ans.diagIm[mRow] = -diagIm[mRow];
1349: return ans;
1350: }
1351:
1352: // TRANSPOSE
1353:
1354: /**
1355: * Returns the transpose of this matrix.
1356: * @return a complex tridiagonal matrix
1357: */
1358: public Matrix transpose() {
1359: final ComplexTridiagonalMatrix ans = new ComplexTridiagonalMatrix(
1360: numRows);
1361: System
1362: .arraycopy(ldiagRe, 1, ans.udiagRe, 0,
1363: ldiagRe.length - 1);
1364: System
1365: .arraycopy(ldiagIm, 1, ans.udiagIm, 0,
1366: ldiagIm.length - 1);
1367: System.arraycopy(diagRe, 0, ans.diagRe, 0, diagRe.length);
1368: System.arraycopy(diagIm, 0, ans.diagIm, 0, diagIm.length);
1369: System
1370: .arraycopy(udiagRe, 0, ans.ldiagRe, 1,
1371: udiagRe.length - 1);
1372: System
1373: .arraycopy(udiagIm, 0, ans.ldiagIm, 1,
1374: udiagIm.length - 1);
1375: return ans;
1376: }
1377:
1378: // MAP ELEMENTS
1379:
1380: /**
1381: * Applies a function on all the matrix elements.
1382: * @param f a user-defined function
1383: * @return a complex tridiagonal matrix
1384: */
1385: public AbstractComplexMatrix mapElements(final ComplexMapping f) {
1386: Complex zeroValue = f.map(Complex.ZERO);
1387: if (zeroValue.mod() <= JSci.GlobalSettings.ZERO_TOL)
1388: return tridiagonalMap(f);
1389: else
1390: return generalMap(f, zeroValue);
1391: }
1392:
1393: private AbstractComplexMatrix tridiagonalMap(ComplexMapping f) {
1394: int mRow = numRows;
1395: final ComplexTridiagonalMatrix ans = new ComplexTridiagonalMatrix(
1396: mRow);
1397: ans.setElement(0, 0, f.map(diagRe[0], diagIm[0]));
1398: ans.setElement(0, 1, f.map(udiagRe[0], udiagIm[0]));
1399: mRow--;
1400: for (int i = 1; i < mRow; i++) {
1401: ans.setElement(i, i - 1, f.map(ldiagRe[i], ldiagIm[i]));
1402: ans.setElement(i, i, f.map(diagRe[i], diagIm[i]));
1403: ans.setElement(i, i + 1, f.map(udiagRe[i], udiagIm[i]));
1404: }
1405: ans.setElement(mRow, mRow - 1, f.map(ldiagRe[mRow],
1406: ldiagIm[mRow]));
1407: ans.setElement(mRow, mRow, f.map(diagRe[mRow], diagIm[mRow]));
1408: return ans;
1409: }
1410:
1411: private AbstractComplexMatrix generalMap(ComplexMapping f,
1412: Complex zeroValue) {
1413: final double arrayRe[][] = new double[numRows][numRows];
1414: final double arrayIm[][] = new double[numRows][numRows];
1415: for (int i = 0; i < numRows; i++) {
1416: for (int j = 0; j < numRows; j++) {
1417: arrayRe[i][j] = zeroValue.real();
1418: arrayIm[i][j] = zeroValue.imag();
1419: }
1420: }
1421: int mRow = numRows;
1422: Complex z = f.map(diagRe[0], diagIm[0]);
1423: arrayRe[0][0] = z.real();
1424: arrayIm[0][0] = z.imag();
1425: z = f.map(udiagRe[0], udiagIm[0]);
1426: arrayRe[0][1] = z.real();
1427: arrayIm[0][1] = z.imag();
1428: mRow--;
1429: for (int i = 1; i < mRow; i++) {
1430: z = f.map(ldiagRe[i], ldiagIm[i]);
1431: arrayRe[i][i - 1] = z.real();
1432: arrayIm[i][i - 1] = z.imag();
1433: z = f.map(diagRe[i], diagIm[i]);
1434: arrayRe[i][i] = z.real();
1435: arrayIm[i][i] = z.imag();
1436: z = f.map(udiagRe[i], udiagIm[i]);
1437: arrayRe[i][i + 1] = z.real();
1438: arrayIm[i][i + 1] = z.imag();
1439: }
1440: z = f.map(ldiagRe[mRow], ldiagIm[mRow]);
1441: arrayRe[mRow][mRow - 1] = z.real();
1442: arrayIm[mRow][mRow - 1] = z.imag();
1443: z = f.map(diagRe[mRow], diagIm[mRow]);
1444: arrayRe[mRow][mRow] = z.real();
1445: arrayIm[mRow][mRow] = z.imag();
1446: return new ComplexSquareMatrix(arrayRe, arrayIm);
1447: }
1448: }
|