Author: mikl
Date: Mon Mar 21 08:14:17 2011
New Revision: 1083698
URL: http://svn.apache.org/viewvc?rev=1083698&view=rev
Log:
Fixes MATH-437
Modified:
commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/AbstractRealMatrix.java
commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/RealMatrix.java
commons/proper/math/trunk/src/test/java/org/apache/commons/math/linear/Array2DRowRealMatrixTest.java
Modified:
commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/AbstractRealMatrix.java
URL:
http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/AbstractRealMatrix.java?rev=1083698&r1=1083697&r2=1083698&view=diff
==============================================================================
---
commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/AbstractRealMatrix.java
(original)
+++
commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/AbstractRealMatrix.java
Mon Mar 21 08:14:17 2011
@@ -17,6 +17,8 @@
package org.apache.commons.math.linear;
+import java.util.ArrayList;
+
import org.apache.commons.math.exception.NoDataException;
import org.apache.commons.math.exception.NotStrictlyPositiveException;
import org.apache.commons.math.exception.DimensionMismatchException;
@@ -147,6 +149,66 @@ public abstract class AbstractRealMatrix
/** {@inheritDoc} */
public RealMatrix preMultiply(final RealMatrix m) {
return m.multiply(this);
+ }
+
+ /** {@inheritDoc} */
+ @Override
+ public RealMatrix power(final int p) {
+ if (p < 0) {
+ throw new IllegalArgumentException("p must be >= 0");
+ }
+
+ if (!isSquare()) {
+ throw new NonSquareMatrixException(getRowDimension(),
getColumnDimension());
+ }
+
+ if (p == 0) {
+ return
MatrixUtils.createRealIdentityMatrix(this.getRowDimension());
+ }
+
+ if (p == 1) {
+ return this.copy();
+ }
+
+ final int power = p - 1;
+
+ /*
+ * Only log_2(p) operations is used by doing as follows:
+ * 5^214 = 5^128 * 5^64 * 5^16 * 5^4 * 5^2
+ *
+ * In general, the same approach is used for A^p.
+ */
+
+ final char[] binaryRepresentation =
Integer.toBinaryString(power).toCharArray();
+ final ArrayList<Integer> nonZeroPositions = new ArrayList<Integer>();
+ int maxI = -1;
+
+ for (int i = 0; i < binaryRepresentation.length; ++i) {
+ if (binaryRepresentation[i] == '1') {
+ final int pos = binaryRepresentation.length - i - 1;
+ nonZeroPositions.add(pos);
+
+ // The positions are taken in turn, so maxI is only changed
once
+ if (maxI == -1) {
+ maxI = pos;
+ }
+ }
+ }
+
+ RealMatrix[] results = new RealMatrix[maxI + 1];
+ results[0] = this.copy();
+
+ for (int i = 1; i <= maxI; ++i) {
+ results[i] = results[i-1].multiply(results[i-1]);
+ }
+
+ RealMatrix result = this.copy();
+
+ for (Integer i : nonZeroPositions) {
+ result = result.multiply(results[i]);
+ }
+
+ return result;
}
/** {@inheritDoc} */
Modified:
commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/RealMatrix.java
URL:
http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/RealMatrix.java?rev=1083698&r1=1083697&r2=1083698&view=diff
==============================================================================
---
commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/RealMatrix.java
(original)
+++
commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/RealMatrix.java
Mon Mar 21 08:14:17 2011
@@ -98,6 +98,16 @@ public interface RealMatrix extends AnyM
* if rowDimension(this) != columnDimension(m)
*/
RealMatrix preMultiply(RealMatrix m);
+
+ /**
+ * Returns the result multiplying this with itself <code>p</code> times.
+ * Depending on the underlying storage, instability for high powers might
occur.
+ * @param p raise this to power p
+ * @return this^p
+ * @throws IllegalArgumentException if p < 0
+ * NonSquareMatrixException if the matrix is not square
+ */
+ RealMatrix power(final int p);
/**
* Returns matrix entries as a two-dimensional array.
Modified:
commons/proper/math/trunk/src/test/java/org/apache/commons/math/linear/Array2DRowRealMatrixTest.java
URL:
http://svn.apache.org/viewvc/commons/proper/math/trunk/src/test/java/org/apache/commons/math/linear/Array2DRowRealMatrixTest.java?rev=1083698&r1=1083697&r2=1083698&view=diff
==============================================================================
---
commons/proper/math/trunk/src/test/java/org/apache/commons/math/linear/Array2DRowRealMatrixTest.java
(original)
+++
commons/proper/math/trunk/src/test/java/org/apache/commons/math/linear/Array2DRowRealMatrixTest.java
Mon Mar 21 08:14:17 2011
@@ -96,6 +96,7 @@ public final class Array2DRowRealMatrixT
// tolerances
protected double entryTolerance = 10E-16;
protected double normTolerance = 10E-14;
+ protected double powerTolerance = 10E-16;
/** test dimensions */
@Test
@@ -222,6 +223,59 @@ public final class Array2DRowRealMatrixT
TestUtils.assertEquals("m3*m4=m5", m3.multiply(m4), m5, entryTolerance);
}
+ public void testPower() {
+ Array2DRowRealMatrix m = new Array2DRowRealMatrix(testData);
+ Array2DRowRealMatrix mInv = new Array2DRowRealMatrix(testDataInv);
+ Array2DRowRealMatrix mPlusInv = new
Array2DRowRealMatrix(testDataPlusInv);
+ Array2DRowRealMatrix identity = new Array2DRowRealMatrix(id);
+
+ TestUtils.assertEquals("m^0", m.power(0),
+ identity, entryTolerance);
+ TestUtils.assertEquals("mInv^0", mInv.power(0),
+ identity, entryTolerance);
+ TestUtils.assertEquals("mPlusInv^0", mPlusInv.power(0),
+ identity, entryTolerance);
+
+ TestUtils.assertEquals("m^1", m.power(1),
+ m, entryTolerance);
+ TestUtils.assertEquals("mInv^1", mInv.power(1),
+ mInv, entryTolerance);
+ TestUtils.assertEquals("mPlusInv^1", mPlusInv.power(1),
+ mPlusInv, entryTolerance);
+
+ RealMatrix C1 = m.copy();
+ RealMatrix C2 = mInv.copy();
+ RealMatrix C3 = mPlusInv.copy();
+
+ for (int i = 2; i <= 10; ++i) {
+ C1 = C1.multiply(m);
+ C2 = C2.multiply(mInv);
+ C3 = C3.multiply(mPlusInv);
+
+ TestUtils.assertEquals("m^" + i, m.power(i),
+ C1, entryTolerance);
+ TestUtils.assertEquals("mInv^" + i, mInv.power(i),
+ C2, entryTolerance);
+ TestUtils.assertEquals("mPlusInv^" + i, mPlusInv.power(i),
+ C3, entryTolerance);
+ }
+
+ try {
+ Array2DRowRealMatrix mNotSquare = new
Array2DRowRealMatrix(testData2T);
+ mNotSquare.power(2);
+ Assert.fail("Expecting NonSquareMatrixException");
+ } catch (NonSquareMatrixException ex) {
+ // ignored
+ }
+
+ try {
+ m.power(-1);
+ Assert.fail("Expecting IllegalArgumentException");
+ } catch (IllegalArgumentException ex) {
+ // ignored
+ }
+ }
+
/** test trace */
@Test
public void testTrace() {
