001: /*
002: * JScience - Java(TM) Tools and Libraries for the Advancement of Sciences.
003: * Copyright (C) 2006 - JScience (http://jscience.org/)
004: * All rights reserved.
005: *
006: * Permission to use, copy, modify, and distribute this software is
007: * freely granted, provided that this notice is preserved.
008: */
009: package org.jscience.mathematics.vector;
010:
011: import java.util.Comparator;
012:
013: import javolution.context.StackContext;
014: import javolution.lang.MathLib;
015: import javolution.lang.Realtime;
016: import javolution.lang.ValueType;
017: import javolution.text.Text;
018: import javolution.text.TextBuilder;
019: import javolution.xml.XMLFormat;
020: import javolution.xml.stream.XMLStreamException;
021:
022: import org.jscience.mathematics.structure.Field;
023: import org.jscience.mathematics.structure.Ring;
024: import org.jscience.mathematics.structure.VectorSpace;
025:
026: /**
027: * <p> This class represents a rectangular table of elements of a ring-like
028: * algebraic structure.</p>
029: *
030: * <p> Instances of this class can be used to resolve system of linear equations
031: * involving <i>any kind</i> of {@link Field Field} elements
032: * (e.g. {@link org.jscience.mathematics.number.Real Real},
033: * {@link org.jscience.mathematics.number.Complex Complex},
034: * {@link org.jscience.physics.amount.Amount Amount<?>},
035: * {@link org.jscience.mathematics.function.Function Function}, etc).
036: * For example:[code]
037: * // Creates a dense matrix (2x2) of Rational numbers.
038: * DenseMatrix<Rational> M = DenseMatrix.valueOf(
039: * { Rational.valueOf(23, 45), Rational.valueOf(33, 75) },
040: * { Rational.valueOf(15, 31), Rational.valueOf(-20, 45)});
041: *
042: * // Creates a sparse matrix (16x2) of Real numbers.
043: * SparseMatrix<Real> M = SparseMatrix.valueOf(
044: * SparseVector.valueOf(16, Real.ZERO, 0, Real.valueOf(5)),
045: * SparseVector.valueOf(16, Real.ZERO, 15, Real.valueOf(-3)));
046: *
047: * // Creates a floating-point (64 bits) matrix (3x2).
048: * Float64Matrix M = Float64Matrix.valueOf(
049: * {{ 1.0, 2.0, 3.0}, { 4.0, 5.0, 6.0}});
050: *
051: * // Creates a complex single column matrix (1x2).
052: * ComplexMatrix M = ComplexMatrix.valueOf(
053: * {{ Complex.valueOf(1.0, 2.0), Complex.valueOf(4.0, 5.0)}}).transpose();
054: *
055: * // Creates an identity matrix (2x2) for modulo integer.
056: * SparseMatrix<ModuloInteger> IDENTITY = SparseMatrix.valueOf(
057: * DenseVector.valueOf(ModuloInteger.ONE, ModuloInteger.ONE), ModuloInteger.ZERO);
058: * [/code]</p>
059: *
060: * <p> Non-commutative field multiplication is supported. Invertible square
061: * matrices may form a non-commutative field (also called a division
062: * ring). In which case this class may be used to resolve system of linear
063: * equations with matrix coefficients.</p>
064: *
065: * <p> Implementation Note: Matrices may use {@link
066: * javolution.context.StackContext StackContext} and {@link
067: * javolution.context.ConcurrentContext ConcurrentContext} in order to
068: * minimize heap allocation and accelerate calculations on multi-core
069: * systems.</p>
070: *
071: * @author <a href="mailto:jean-marie@dautelle.com">Jean-Marie Dautelle</a>
072: * @version 3.3, December 24, 2006
073: * @see <a href="http://en.wikipedia.org/wiki/Matrix_%28mathematics%29">
074: * Wikipedia: Matrix (mathematics)</a>
075: */
076: public abstract class Matrix<F extends Field<F>> implements
077: VectorSpace<Matrix<F>, F>, Ring<Matrix<F>>, ValueType, Realtime {
078:
079: /**
080: * Holds the default XML representation for matrices. For example:[code]
081: * <DenseMatrix rows="2" columns="2">
082: * <Complex real="1.0" imaginary="0.0" />
083: * <Complex real="0.0" imaginary="1.0" />
084: * <Complex real="0.0" imaginary="0.4" />
085: * <Complex real="-5.0" imaginary="-1.0" />
086: * </DenseMatrix>[/code]
087: */
088: @SuppressWarnings("unchecked")
089: protected static final XMLFormat<Matrix> XML = new XMLFormat<Matrix>(
090: Matrix.class) {
091:
092: @Override
093: public void read(InputElement xml, Matrix M)
094: throws XMLStreamException {
095: // Nothing to do.
096: }
097:
098: @Override
099: public void write(Matrix M, OutputElement xml)
100: throws XMLStreamException {
101: final int m = M.getNumberOfRows();
102: final int n = M.getNumberOfColumns();
103: xml.setAttribute("rows", m);
104: xml.setAttribute("columns", n);
105: for (int i = 0; i < m; i++) {
106: for (int j = 0; j < n; j++) {
107: xml.add(M.get(i, j));
108: }
109: }
110: }
111: };
112:
113: /**
114: * Default constructor (for sub-classes).
115: */
116: protected Matrix() {
117: }
118:
119: /**
120: * Returns the number of rows <code>m</code> for this matrix.
121: *
122: * @return m, the number of rows.
123: */
124: public abstract int getNumberOfRows();
125:
126: /**
127: * Returns the number of columns <code>n</code> for this matrix.
128: *
129: * @return n, the number of columns.
130: */
131: public abstract int getNumberOfColumns();
132:
133: /**
134: * Returns a single element from this matrix.
135: *
136: * @param i the row index (range [0..m[).
137: * @param j the column index (range [0..n[).
138: * @return the element read at [i,j].
139: * @throws IndexOutOfBoundsException <code>
140: * ((i < 0) || (i >= m)) || ((j < 0) || (j >= n))</code>
141: */
142: public abstract F get(int i, int j);
143:
144: /**
145: * Returns the row identified by the specified index in this matrix.
146: *
147: * @param i the row index (range [0..m[).
148: * @return the vector holding the specified row.
149: * @throws IndexOutOfBoundsException <code>(i < 0) || (i >= m)</code>
150: */
151: public abstract Vector<F> getRow(int i);
152:
153: /**
154: * Returns the column identified by the specified index in this matrix.
155: *
156: * @param j the column index (range [0..n[).
157: * @return the vector holding the specified column.
158: * @throws IndexOutOfBoundsException <code>(j < 0) || (j >= n)</code>
159: */
160: public abstract Vector<F> getColumn(int j);
161:
162: /**
163: * Returns the diagonal vector.
164: *
165: * @return the vector holding the diagonal elements.
166: */
167: public abstract Vector<F> getDiagonal();
168:
169: /**
170: * Returns the negation of this matrix.
171: *
172: * @return <code>-this</code>.
173: */
174: public abstract Matrix<F> opposite();
175:
176: /**
177: * Returns the sum of this matrix with the one specified.
178: *
179: * @param that the matrix to be added.
180: * @return <code>this + that</code>.
181: * @throws DimensionException matrices's dimensions are different.
182: */
183: public abstract Matrix<F> plus(Matrix<F> that);
184:
185: /**
186: * Returns the difference between this matrix and the one specified.
187: *
188: * @param that the matrix to be subtracted.
189: * @return <code>this - that</code>.
190: * @throws DimensionException matrices's dimensions are different.
191: */
192: public Matrix<F> minus(Matrix<F> that) {
193: return this .plus(that.opposite());
194: }
195:
196: /**
197: * Returns the product of this matrix by the specified factor.
198: *
199: * @param k the coefficient multiplier.
200: * @return <code>this · k</code>
201: */
202: public abstract Matrix<F> times(F k);
203:
204: /**
205: * Returns the product of this matrix by the specified vector.
206: *
207: * @param v the vector.
208: * @return <code>this · v</code>
209: * @throws DimensionException if <code>
210: * v.getDimension() != this.getNumberOfColumns()<code>
211: */
212: public abstract Vector<F> times(Vector<F> v);
213:
214: /**
215: * Returns the product of this matrix with the one specified.
216: *
217: * @param that the matrix multiplier.
218: * @return <code>this · that</code>.
219: * @throws DimensionException if <code>
220: * this.getNumberOfColumns() != that.getNumberOfRows()</code>.
221: */
222: public abstract Matrix<F> times(Matrix<F> that);
223:
224: /**
225: * Returns the inverse of this matrix (must be square).
226: *
227: * @return <code>1 / this</code>
228: * @throws DimensionException if this matrix is not square.
229: */
230: public abstract Matrix<F> inverse();
231:
232: /**
233: * Returns this matrix divided by the one specified.
234: *
235: * @param that the matrix divisor.
236: * @return <code>this / that</code>.
237: * @throws DimensionException if that matrix is not square or dimensions
238: * do not match.
239: */
240: public Matrix<F> divide(Matrix<F> that) {
241: return this .times(that.inverse());
242: }
243:
244: /**
245: * Returns the inverse or pseudo-inverse if this matrix if not square.
246: *
247: * <p> Note: To resolve the equation <code>A * X = B</code>,
248: * it is usually faster to calculate <code>A.lu().solve(B)</code>
249: * rather than <code>A.inverse().times(B)</code>.</p>
250: *
251: * @return the inverse or pseudo-inverse of this matrix.
252: */
253: public Matrix<F> pseudoInverse() {
254: if (isSquare())
255: return this .inverse();
256: Matrix<F> this Transpose = this .transpose();
257: return (this Transpose.times(this )).inverse().times(
258: this Transpose);
259: }
260:
261: /**
262: * Returns the determinant of this matrix.
263: *
264: * @return this matrix determinant.
265: * @throws DimensionException if this matrix is not square.
266: */
267: public abstract F determinant();
268:
269: /**
270: * Returns the transpose of this matrix.
271: *
272: * @return <code>A'</code>.
273: */
274: public abstract Matrix<F> transpose();
275:
276: /**
277: * Returns the cofactor of an element in this matrix. It is the value
278: * obtained by evaluating the determinant formed by the elements not in
279: * that particular row or column.
280: *
281: * @param i the row index.
282: * @param j the column index.
283: * @return the cofactor of <code>THIS[i,j]</code>.
284: * @throws DimensionException matrix is not square or its dimension
285: * is less than 2.
286: */
287: public abstract F cofactor(int i, int j);
288:
289: /**
290: * Returns the adjoint of this matrix. It is obtained by replacing each
291: * element in this matrix with its cofactor and applying a + or - sign
292: * according (-1)**(i+j), and then finding the transpose of the resulting
293: * matrix.
294: *
295: * @return the adjoint of this matrix.
296: * @throws DimensionException if this matrix is not square or if
297: * its dimension is less than 2.
298: */
299: public abstract Matrix<F> adjoint();
300:
301: /**
302: * Indicates if this matrix is square.
303: *
304: * @return <code>getNumberOfRows() == getNumberOfColumns()</code>
305: */
306: public boolean isSquare() {
307: return getNumberOfRows() == getNumberOfColumns();
308: }
309:
310: /**
311: * Solves this matrix for the specified vector (returns <code>x</code>
312: * such as <code>this · x = y</code>).
313: *
314: * @param y the vector for which the solution is calculated.
315: * @return <code>x</code> such as <code>this · x = y</code>
316: * @throws DimensionException if that matrix is not square or dimensions
317: * do not match.
318: */
319: public Vector<F> solve(Vector<F> y) {
320: DenseMatrix<F> M = DenseMatrix.newInstance(y.getDimension(),
321: true);
322: M._rows.add(DenseVector.valueOf(y));
323: return solve(M).getColumn(0);
324: }
325:
326: /**
327: * Solves this matrix for the specified matrix (returns <code>x</code>
328: * such as <code>this · x = y</code>).
329: *
330: * @param y the matrix for which the solution is calculated.
331: * @return <code>x</code> such as <code>this · x = y</code>
332: * @throws DimensionException if that matrix is not square or dimensions
333: * do not match.
334: */
335: public Matrix<F> solve(Matrix<F> y) {
336: return LUDecomposition.valueOf(this ).solve(y); // Default implementation.
337: }
338:
339: /**
340: * Returns this matrix raised at the specified exponent.
341: *
342: * @param exp the exponent.
343: * @return <code>this<sup>exp</sup></code>
344: * @throws DimensionException if this matrix is not square.
345: */
346: public Matrix<F> pow(int exp) {
347: if (exp > 0) {
348: StackContext.enter();
349: try {
350: Matrix<F> pow2 = this ;
351: Matrix<F> result = null;
352: while (exp >= 1) { // Iteration.
353: if ((exp & 1) == 1) {
354: result = (result == null) ? pow2 : result
355: .times(pow2);
356: }
357: pow2 = pow2.times(pow2);
358: exp >>>= 1;
359: }
360: return StackContext.outerCopy(result);
361: } finally {
362: StackContext.exit();
363: }
364: } else if (exp == 0) {
365: return this .times(this .inverse()); // Identity.
366: } else {
367: return this .pow(-exp).inverse();
368: }
369: }
370:
371: /**
372: * Returns the trace of this matrix.
373: *
374: * @return the sum of the diagonal elements.
375: */
376: public F trace() {
377: F sum = this .get(0, 0);
378: for (int i = MathLib.min(getNumberOfColumns(),
379: getNumberOfRows()); --i > 0;) {
380: sum = sum.plus(get(i, i));
381: }
382: return sum;
383: }
384:
385: /**
386: * Returns the linear algebraic matrix tensor product of this matrix
387: * and another (Kronecker product). The default implementation returns
388: * a {@link DenseMatrix}.
389: *
390: * @param that the second matrix.
391: * @return <code>this ⊗ that</code>
392: * @see <a href="http://en.wikipedia.org/wiki/Kronecker_product">
393: * Wikipedia: Kronecker Product</a>
394: */
395: public abstract Matrix<F> tensor(Matrix<F> that);
396:
397: /**
398: * Returns the vectorization of this matrix. The vectorization of
399: * a matrix is the column vector obtain by stacking the columns of the
400: * matrix on top of one another. The default implementation returns
401: * a {@link DenseVector}.
402: *
403: * @return the vectorization of this matrix.
404: * @see <a href="http://en.wikipedia.org/wiki/Vectorization_%28mathematics%29">
405: * Wikipedia: Vectorization.</a>
406: */
407: public abstract Vector<F> vectorization();
408:
409: /**
410: * Returns the text representation of this matrix.
411: *
412: * @return the text representation of this matrix.
413: */
414: public Text toText() {
415: final int m = this .getNumberOfRows();
416: final int n = this .getNumberOfColumns();
417: TextBuilder tmp = TextBuilder.newInstance();
418: tmp.append('{');
419: for (int i = 0; i < m; i++) {
420: tmp.append('{');
421: for (int j = 0; j < n; j++) {
422: tmp.append(get(i, j));
423: if (j != n - 1) {
424: tmp.append(", ");
425: }
426: }
427: tmp.append("}");
428: if (i != m - 1) {
429: tmp.append(",\n");
430: }
431: }
432: tmp.append("}");
433: Text txt = tmp.toText();
434: TextBuilder.recycle(tmp);
435: return txt;
436: }
437:
438: /**
439: * Returns the text representation of this matrix as a
440: * <code>java.lang.String</code>.
441: *
442: * @return <code>toText().toString()</code>
443: */
444: public final String toString() {
445: return toText().toString();
446: }
447:
448: /**
449: * Indicates if this matrix can be considered equals to the one
450: * specified using the specified comparator when testing for
451: * element equality. The specified comparator may allow for some
452: * tolerance in the difference between the matrix elements.
453: *
454: * @param that the matrix to compare for equality.
455: * @param cmp the comparator to use when testing for element equality.
456: * @return <code>true</code> if this matrix and the specified matrix are
457: * both matrices with equal elements according to the specified
458: * comparator; <code>false</code> otherwise.
459: */
460: public boolean equals(Matrix<F> that, Comparator<F> cmp) {
461: if (this == that)
462: return true;
463: final int m = this .getNumberOfRows();
464: final int n = this .getNumberOfColumns();
465: if ((that.getNumberOfRows() != m)
466: || (that.getNumberOfColumns() != n))
467: return false;
468: for (int i = m; --i >= 0;) {
469: for (int j = n; --j >= 0;) {
470: if (cmp.compare(this .get(i, j), that.get(i, j)) != 0)
471: return false;
472: }
473: }
474: return true;
475: }
476:
477: /**
478: * Indicates if this matrix is strictly equal to the object specified.
479: *
480: * @param that the object to compare for equality.
481: * @return <code>true</code> if this matrix and the specified object are
482: * both matrices with equal elements; <code>false</code> otherwise.
483: * @see #equals(Matrix, Comparator)
484: */
485: public boolean equals(Object that) {
486: if (this == that)
487: return true;
488: if (!(that instanceof Matrix))
489: return false;
490: final int m = this .getNumberOfRows();
491: final int n = this .getNumberOfColumns();
492: Matrix<?> M = (Matrix<?>) that;
493: if ((M.getNumberOfRows() != m) || (M.getNumberOfColumns() != n))
494: return false;
495: for (int i = m; --i >= 0;) {
496: for (int j = n; --j >= 0;) {
497: if (!this .get(i, j).equals(M.get(i, j)))
498: return false;
499: }
500: }
501: return true;
502: }
503:
504: /**
505: * Returns a hash code value for this matrix.
506: * Equals objects have equal hash codes.
507: *
508: * @return this matrix hash code value.
509: * @see #equals
510: */
511: public int hashCode() {
512: final int m = this .getNumberOfRows();
513: final int n = this .getNumberOfColumns();
514: int code = 0;
515: for (int i = m; --i >= 0;) {
516: for (int j = n; --j >= 0;) {
517: code += get(i, j).hashCode();
518: }
519: }
520: return code;
521: }
522:
523: /**
524: * Returns a copy of this matrix
525: * {@link javolution.context.AllocatorContext allocated}
526: * by the calling thread (possibly on the stack).
527: *
528: * @return an identical and independant copy of this matrix.
529: */
530: public abstract Matrix<F> copy();
531:
532: }
|