/*
 * Licensed to the Apache Software Foundation (ASF) under one or more
 * contributor license agreements.  See the NOTICE file distributed with
 * this work for additional information regarding copyright ownership.
 * The ASF licenses this file to You under the Apache License, Version 2.0
 * (the "License"); you may not use this file except in compliance with
 * the License.  You may obtain a copy of the License at
 *
 *      http://www.apache.org/licenses/LICENSE-2.0
 *
 * Unless required by applicable law or agreed to in writing, software
 * distributed under the License is distributed on an "AS IS" BASIS,
 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 * See the License for the specific language governing permissions and
 * limitations under the License.
 */

package org.apache.commons.math.linear;

import java.io.Serializable;
import java.math.BigDecimal;

import org.apache.commons.math.MathRuntimeException;

/**
 * Implementation of {@link BigMatrix} using a BigDecimal[][] array to store
 * entries and
 * <a href="http://www.math.gatech.edu/~bourbaki/math2601/Web-notes/2num.pdf">
 * LU decompostion</a> to support linear system solution and inverse.
 * <p>
 * The LU decompostion is performed as needed, to support the following
 * operations:
 * <ul>
 * <li>solve</li>
 * <li>isSingular</li>
 * <li>getDeterminant</li>
 * <li>inverse</li>
 * </ul>
 * </p>
 * <p>
 * <strong>Usage notes</strong>:<br>
 * <ul>
 * <li>The LU decomposition is stored and reused on subsequent calls. If matrix
 * data are modified using any of the public setXxx methods, the saved
 * decomposition is discarded. If data are modified via references to the
 * underlying array obtained using <code>getDataRef()</code>, then the stored LU
 * decomposition will not be discarded. In this case, you need to explicitly
 * invoke <code>LUDecompose()</code> to recompute the decomposition before using
 * any of the methods above.</li>
 * <li>As specified in the {@link BigMatrix} interface, matrix element indexing
 * is 0-based -- e.g., <code>getEntry(0, 0)</code> returns the element in the
 * first row, first column of the matrix.</li>
 * </ul>
 * </p>
 * 
 * @deprecated as of 2.0, replaced by {@link Array2DRowFieldMatrix} with a
 *             {@link org.apache.commons.math.util.BigReal} parameter
 * @version $Revision$ $Date$
 */
@Deprecated
public class BigMatrixImpl implements BigMatrix, Serializable {

	/** Serialization id */
	private static final long serialVersionUID = -1011428905656140431L;

	/** Entries of the matrix */
	protected BigDecimal data[][] = null;

	/**
	 * Entries of cached LU decomposition. All updates to data (other than
	 * luDecompose()) *must* set this to null
	 */
	protected BigDecimal lu[][] = null;

	/** Permutation associated with LU decomposition */
	protected int[] permutation = null;

	/** Parity of the permutation associated with the LU decomposition */
	protected int parity = 1;

	/** Rounding mode for divisions **/
	private int roundingMode = BigDecimal.ROUND_HALF_UP;

	/*** BigDecimal scale ***/
	private int scale = 64;

	/** Bound to determine effective singularity in LU decomposition */
	private static final BigDecimal TOO_SMALL = new BigDecimal(10E-12);

	/** BigDecimal 0 */
	static final BigDecimal ZERO = new BigDecimal(0);
	/** BigDecimal 1 */
	static final BigDecimal ONE = new BigDecimal(1);

	/**
	 * Creates a matrix with no data
	 */
	public BigMatrixImpl() {
	}

	/**
	 * Create a new BigMatrix with the supplied row and column dimensions.
	 *
	 * @param rowDimension
	 *            the number of rows in the new matrix
	 * @param columnDimension
	 *            the number of columns in the new matrix
	 * @throws IllegalArgumentException
	 *             if row or column dimension is not positive
	 */
	public BigMatrixImpl(int rowDimension, int columnDimension) {
		if (rowDimension <= 0) {
			throw MathRuntimeException.createIllegalArgumentException("invalid row dimension {0} (must be positive)",
					rowDimension);
		}
		if (columnDimension <= 0) {
			throw MathRuntimeException.createIllegalArgumentException("invalid column dimension {0} (must be positive)",
					columnDimension);
		}
		data = new BigDecimal[rowDimension][columnDimension];
		lu = null;
	}

	/**
	 * Create a new BigMatrix using <code>d</code> as the underlying data array.
	 * <p>
	 * The input array is copied, not referenced. This constructor has the same
	 * effect as calling {@link #BigMatrixImpl(BigDecimal[][], boolean)} with
	 * the second argument set to <code>true</code>.
	 * </p>
	 *
	 * @param d
	 *            data for new matrix
	 * @throws IllegalArgumentException
	 *             if <code>d</code> is not rectangular (not all rows have the
	 *             same length) or empty
	 * @throws NullPointerException
	 *             if <code>d</code> is null
	 */
	public BigMatrixImpl(BigDecimal[][] d) {
		this.copyIn(d);
		lu = null;
	}

	/**
	 * Create a new BigMatrix using the input array as the underlying data
	 * array.
	 * <p>
	 * If an array is built specially in order to be embedded in a BigMatrix and
	 * not used directly, the <code>copyArray</code> may be set to
	 * <code>false</code. This will prevent the copying and improve
	 * performance as no new array will be built and no data will be copied.</p>
	 * &#64;param d data for new matrix
	 * &#64;param copyArray if true, the input array will be copied, otherwise
	 * it will be referenced
	 * @throws IllegalArgumentException if <code>d</code> is not rectangular
	 *             (not all rows have the same length) or empty
	 * 
	 * @throws NullPointerException
	 *             if <code>d</code> is null
	 * @see #BigMatrixImpl(BigDecimal[][])
	 */
	public BigMatrixImpl(BigDecimal[][] d, boolean copyArray) {
		if (copyArray) {
			copyIn(d);
		} else {
			if (d == null) {
				throw new NullPointerException();
			}
			final int nRows = d.length;
			if (nRows == 0) {
				throw MathRuntimeException.createIllegalArgumentException("matrix must have at least one row");
			}

			final int nCols = d[0].length;
			if (nCols == 0) {
				throw MathRuntimeException.createIllegalArgumentException("matrix must have at least one column");
			}
			for (int r = 1; r < nRows; r++) {
				if (d[r].length != nCols) {
					throw MathRuntimeException.createIllegalArgumentException(
							"some rows have length {0} while others have length {1}", nCols, d[r].length);
				}
			}
			data = d;
		}
		lu = null;
	}

	/**
	 * Create a new BigMatrix using <code>d</code> as the underlying data array.
	 * <p>
	 * Since the underlying array will hold <code>BigDecimal</code> instances,
	 * it will be created.
	 * </p>
	 *
	 * @param d
	 *            data for new matrix
	 * @throws IllegalArgumentException
	 *             if <code>d</code> is not rectangular (not all rows have the
	 *             same length) or empty
	 * @throws NullPointerException
	 *             if <code>d</code> is null
	 */
	public BigMatrixImpl(double[][] d) {
		final int nRows = d.length;
		if (nRows == 0) {
			throw MathRuntimeException.createIllegalArgumentException("matrix must have at least one row");
		}

		final int nCols = d[0].length;
		if (nCols == 0) {
			throw MathRuntimeException.createIllegalArgumentException("matrix must have at least one column");
		}
		for (int row = 1; row < nRows; row++) {
			if (d[row].length != nCols) {
				throw MathRuntimeException.createIllegalArgumentException(
						"some rows have length {0} while others have length {1}", nCols, d[row].length);
			}
		}
		this.copyIn(d);
		lu = null;
	}

	/**
	 * Create a new BigMatrix using the values represented by the strings in
	 * <code>d</code> as the underlying data array.
	 *
	 * @param d
	 *            data for new matrix
	 * @throws IllegalArgumentException
	 *             if <code>d</code> is not rectangular (not all rows have the
	 *             same length) or empty
	 * @throws NullPointerException
	 *             if <code>d</code> is null
	 */
	public BigMatrixImpl(String[][] d) {
		final int nRows = d.length;
		if (nRows == 0) {
			throw MathRuntimeException.createIllegalArgumentException("matrix must have at least one row");
		}

		final int nCols = d[0].length;
		if (nCols == 0) {
			throw MathRuntimeException.createIllegalArgumentException("matrix must have at least one column");
		}
		for (int row = 1; row < nRows; row++) {
			if (d[row].length != nCols) {
				throw MathRuntimeException.createIllegalArgumentException(
						"some rows have length {0} while others have length {1}", nCols, d[row].length);
			}
		}
		this.copyIn(d);
		lu = null;
	}

	/**
	 * Create a new (column) BigMatrix using <code>v</code> as the data for the
	 * unique column of the <code>v.length x 1</code> matrix created.
	 * <p>
	 * The input array is copied, not referenced.
	 * </p>
	 *
	 * @param v
	 *            column vector holding data for new matrix
	 */
	public BigMatrixImpl(BigDecimal[] v) {
		final int nRows = v.length;
		data = new BigDecimal[nRows][1];
		for (int row = 0; row < nRows; row++) {
			data[row][0] = v[row];
		}
	}

	/**
	 * Create a new BigMatrix which is a copy of this.
	 *
	 * @return the cloned matrix
	 */
	public BigMatrix copy() {
		return new BigMatrixImpl(this.copyOut(), false);
	}

	/**
	 * Compute the sum of this and <code>m</code>.
	 *
	 * @param m
	 *            matrix to be added
	 * @return this + m
	 * @throws IllegalArgumentException
	 *             if m is not the same size as this
	 */
	public BigMatrix add(BigMatrix m) throws IllegalArgumentException {
		try {
			return add((BigMatrixImpl) m);
		} catch (ClassCastException cce) {

			// safety check
			MatrixUtils.checkAdditionCompatible(this, m);

			final int rowCount = getRowDimension();
			final int columnCount = getColumnDimension();
			final BigDecimal[][] outData = new BigDecimal[rowCount][columnCount];
			for (int row = 0; row < rowCount; row++) {
				final BigDecimal[] dataRow = data[row];
				final BigDecimal[] outDataRow = outData[row];
				for (int col = 0; col < columnCount; col++) {
					outDataRow[col] = dataRow[col].add(m.getEntry(row, col));
				}
			}
			return new BigMatrixImpl(outData, false);
		}
	}

	/**
	 * Compute the sum of this and <code>m</code>.
	 *
	 * @param m
	 *            matrix to be added
	 * @return this + m
	 * @throws IllegalArgumentException
	 *             if m is not the same size as this
	 */
	public BigMatrixImpl add(BigMatrixImpl m) throws IllegalArgumentException {

		// safety check
		MatrixUtils.checkAdditionCompatible(this, m);

		final int rowCount = getRowDimension();
		final int columnCount = getColumnDimension();
		final BigDecimal[][] outData = new BigDecimal[rowCount][columnCount];
		for (int row = 0; row < rowCount; row++) {
			final BigDecimal[] dataRow = data[row];
			final BigDecimal[] mRow = m.data[row];
			final BigDecimal[] outDataRow = outData[row];
			for (int col = 0; col < columnCount; col++) {
				outDataRow[col] = dataRow[col].add(mRow[col]);
			}
		}
		return new BigMatrixImpl(outData, false);
	}

	/**
	 * Compute this minus <code>m</code>.
	 *
	 * @param m
	 *            matrix to be subtracted
	 * @return this + m
	 * @throws IllegalArgumentException
	 *             if m is not the same size as this
	 */
	public BigMatrix subtract(BigMatrix m) throws IllegalArgumentException {
		try {
			return subtract((BigMatrixImpl) m);
		} catch (ClassCastException cce) {

			// safety check
			MatrixUtils.checkSubtractionCompatible(this, m);

			final int rowCount = getRowDimension();
			final int columnCount = getColumnDimension();
			final BigDecimal[][] outData = new BigDecimal[rowCount][columnCount];
			for (int row = 0; row < rowCount; row++) {
				final BigDecimal[] dataRow = data[row];
				final BigDecimal[] outDataRow = outData[row];
				for (int col = 0; col < columnCount; col++) {
					outDataRow[col] = dataRow[col].subtract(getEntry(row, col));
				}
			}
			return new BigMatrixImpl(outData, false);
		}
	}

	/**
	 * Compute this minus <code>m</code>.
	 *
	 * @param m
	 *            matrix to be subtracted
	 * @return this + m
	 * @throws IllegalArgumentException
	 *             if m is not the same size as this
	 */
	public BigMatrixImpl subtract(BigMatrixImpl m) throws IllegalArgumentException {

		// safety check
		MatrixUtils.checkSubtractionCompatible(this, m);

		final int rowCount = getRowDimension();
		final int columnCount = getColumnDimension();
		final BigDecimal[][] outData = new BigDecimal[rowCount][columnCount];
		for (int row = 0; row < rowCount; row++) {
			final BigDecimal[] dataRow = data[row];
			final BigDecimal[] mRow = m.data[row];
			final BigDecimal[] outDataRow = outData[row];
			for (int col = 0; col < columnCount; col++) {
				outDataRow[col] = dataRow[col].subtract(mRow[col]);
			}
		}
		return new BigMatrixImpl(outData, false);
	}

	/**
	 * Returns the result of adding d to each entry of this.
	 *
	 * @param d
	 *            value to be added to each entry
	 * @return d + this
	 */
	public BigMatrix scalarAdd(BigDecimal d) {
		final int rowCount = getRowDimension();
		final int columnCount = getColumnDimension();
		final BigDecimal[][] outData = new BigDecimal[rowCount][columnCount];
		for (int row = 0; row < rowCount; row++) {
			final BigDecimal[] dataRow = data[row];
			final BigDecimal[] outDataRow = outData[row];
			for (int col = 0; col < columnCount; col++) {
				outDataRow[col] = dataRow[col].add(d);
			}
		}
		return new BigMatrixImpl(outData, false);
	}

	/**
	 * Returns the result of multiplying each entry of this by <code>d</code>
	 * 
	 * @param d
	 *            value to multiply all entries by
	 * @return d * this
	 */
	public BigMatrix scalarMultiply(BigDecimal d) {
		final int rowCount = getRowDimension();
		final int columnCount = getColumnDimension();
		final BigDecimal[][] outData = new BigDecimal[rowCount][columnCount];
		for (int row = 0; row < rowCount; row++) {
			final BigDecimal[] dataRow = data[row];
			final BigDecimal[] outDataRow = outData[row];
			for (int col = 0; col < columnCount; col++) {
				outDataRow[col] = dataRow[col].multiply(d);
			}
		}
		return new BigMatrixImpl(outData, false);
	}

	/**
	 * Returns the result of postmultiplying this by <code>m</code>.
	 * 
	 * @param m
	 *            matrix to postmultiply by
	 * @return this*m
	 * @throws IllegalArgumentException
	 *             if columnDimension(this) != rowDimension(m)
	 */
	public BigMatrix multiply(BigMatrix m) throws IllegalArgumentException {
		try {
			return multiply((BigMatrixImpl) m);
		} catch (ClassCastException cce) {

			// safety check
			MatrixUtils.checkMultiplicationCompatible(this, m);

			final int nRows = this.getRowDimension();
			final int nCols = m.getColumnDimension();
			final int nSum = this.getColumnDimension();
			final BigDecimal[][] outData = new BigDecimal[nRows][nCols];
			for (int row = 0; row < nRows; row++) {
				final BigDecimal[] dataRow = data[row];
				final BigDecimal[] outDataRow = outData[row];
				for (int col = 0; col < nCols; col++) {
					BigDecimal sum = ZERO;
					for (int i = 0; i < nSum; i++) {
						sum = sum.add(dataRow[i].multiply(m.getEntry(i, col)));
					}
					outDataRow[col] = sum;
				}
			}
			return new BigMatrixImpl(outData, false);
		}
	}

	/**
	 * Returns the result of postmultiplying this by <code>m</code>.
	 * 
	 * @param m
	 *            matrix to postmultiply by
	 * @return this*m
	 * @throws IllegalArgumentException
	 *             if columnDimension(this) != rowDimension(m)
	 */
	public BigMatrixImpl multiply(BigMatrixImpl m) throws IllegalArgumentException {

		// safety check
		MatrixUtils.checkMultiplicationCompatible(this, m);

		final int nRows = this.getRowDimension();
		final int nCols = m.getColumnDimension();
		final int nSum = this.getColumnDimension();
		final BigDecimal[][] outData = new BigDecimal[nRows][nCols];
		for (int row = 0; row < nRows; row++) {
			final BigDecimal[] dataRow = data[row];
			final BigDecimal[] outDataRow = outData[row];
			for (int col = 0; col < nCols; col++) {
				BigDecimal sum = ZERO;
				for (int i = 0; i < nSum; i++) {
					sum = sum.add(dataRow[i].multiply(m.data[i][col]));
				}
				outDataRow[col] = sum;
			}
		}
		return new BigMatrixImpl(outData, false);
	}

	/**
	 * Returns the result premultiplying this by <code>m</code>.
	 * 
	 * @param m
	 *            matrix to premultiply by
	 * @return m * this
	 * @throws IllegalArgumentException
	 *             if rowDimension(this) != columnDimension(m)
	 */
	public BigMatrix preMultiply(BigMatrix m) throws IllegalArgumentException {
		return m.multiply(this);
	}

	/**
	 * Returns matrix entries as a two-dimensional array.
	 * <p>
	 * Makes a fresh copy of the underlying data.
	 * </p>
	 *
	 * @return 2-dimensional array of entries
	 */
	public BigDecimal[][] getData() {
		return copyOut();
	}

	/**
	 * Returns matrix entries as a two-dimensional array.
	 * <p>
	 * Makes a fresh copy of the underlying data converted to
	 * <code>double</code> values.
	 * </p>
	 *
	 * @return 2-dimensional array of entries
	 */
	public double[][] getDataAsDoubleArray() {
		final int nRows = getRowDimension();
		final int nCols = getColumnDimension();
		final double d[][] = new double[nRows][nCols];
		for (int i = 0; i < nRows; i++) {
			for (int j = 0; j < nCols; j++) {
				d[i][j] = data[i][j].doubleValue();
			}
		}
		return d;
	}

	/**
	 * Returns a reference to the underlying data array.
	 * <p>
	 * Does not make a fresh copy of the underlying data.
	 * </p>
	 *
	 * @return 2-dimensional array of entries
	 */
	public BigDecimal[][] getDataRef() {
		return data;
	}

	/***
	 * Gets the rounding mode for division operations The default is
	 * {@link java.math.BigDecimal#ROUND_HALF_UP}
	 * 
	 * @see BigDecimal
	 * @return the rounding mode.
	 */
	public int getRoundingMode() {
		return roundingMode;
	}

	/***
	 * Sets the rounding mode for decimal divisions.
	 * 
	 * @see BigDecimal
	 * @param roundingMode
	 *            rounding mode for decimal divisions
	 */
	public void setRoundingMode(int roundingMode) {
		this.roundingMode = roundingMode;
	}

	/***
	 * Sets the scale for division operations. The default is 64
	 * 
	 * @see BigDecimal
	 * @return the scale
	 */
	public int getScale() {
		return scale;
	}

	/***
	 * Sets the scale for division operations.
	 * 
	 * @see BigDecimal
	 * @param scale
	 *            scale for division operations
	 */
	public void setScale(int scale) {
		this.scale = scale;
	}

	/**
	 * Returns the
	 * <a href="http://mathworld.wolfram.com/MaximumAbsoluteRowSumNorm.html">
	 * maximum absolute row sum norm</a> of the matrix.
	 *
	 * @return norm
	 */
	public BigDecimal getNorm() {
		BigDecimal maxColSum = ZERO;
		for (int col = 0; col < this.getColumnDimension(); col++) {
			BigDecimal sum = ZERO;
			for (int row = 0; row < this.getRowDimension(); row++) {
				sum = sum.add(data[row][col].abs());
			}
			maxColSum = maxColSum.max(sum);
		}
		return maxColSum;
	}

	/**
	 * Gets a submatrix. Rows and columns are indicated counting from 0 to n-1.
	 *
	 * @param startRow
	 *            Initial row index
	 * @param endRow
	 *            Final row index
	 * @param startColumn
	 *            Initial column index
	 * @param endColumn
	 *            Final column index
	 * @return The subMatrix containing the data of the specified rows and
	 *         columns
	 * @exception MatrixIndexException
	 *                if row or column selections are not valid
	 */
	public BigMatrix getSubMatrix(int startRow, int endRow, int startColumn, int endColumn)
			throws MatrixIndexException {

		MatrixUtils.checkRowIndex(this, startRow);
		MatrixUtils.checkRowIndex(this, endRow);
		if (startRow > endRow) {
			throw new MatrixIndexException("initial row {0} after final row {1}", startRow, endRow);
		}

		MatrixUtils.checkColumnIndex(this, startColumn);
		MatrixUtils.checkColumnIndex(this, endColumn);
		if (startColumn > endColumn) {
			throw new MatrixIndexException("initial column {0} after final column {1}", startColumn, endColumn);
		}

		final BigDecimal[][] subMatrixData = new BigDecimal[endRow - startRow + 1][endColumn - startColumn + 1];
		for (int i = startRow; i <= endRow; i++) {
			System.arraycopy(data[i], startColumn, subMatrixData[i - startRow], 0, endColumn - startColumn + 1);
		}

		return new BigMatrixImpl(subMatrixData, false);

	}

	/**
	 * Gets a submatrix. Rows and columns are indicated counting from 0 to n-1.
	 *
	 * @param selectedRows
	 *            Array of row indices must be non-empty
	 * @param selectedColumns
	 *            Array of column indices must be non-empty
	 * @return The subMatrix containing the data in the specified rows and
	 *         columns
	 * @exception MatrixIndexException
	 *                if supplied row or column index arrays are not valid
	 */
	public BigMatrix getSubMatrix(int[] selectedRows, int[] selectedColumns) throws MatrixIndexException {

		if (selectedRows.length * selectedColumns.length == 0) {
			if (selectedRows.length == 0) {
				throw new MatrixIndexException("empty selected row index array");
			}
			throw new MatrixIndexException("empty selected column index array");
		}

		final BigDecimal[][] subMatrixData = new BigDecimal[selectedRows.length][selectedColumns.length];
		try {
			for (int i = 0; i < selectedRows.length; i++) {
				final BigDecimal[] subI = subMatrixData[i];
				final BigDecimal[] dataSelectedI = data[selectedRows[i]];
				for (int j = 0; j < selectedColumns.length; j++) {
					subI[j] = dataSelectedI[selectedColumns[j]];
				}
			}
		} catch (ArrayIndexOutOfBoundsException e) {
			// we redo the loop with checks enabled
			// in order to generate an appropriate message
			for (final int row : selectedRows) {
				MatrixUtils.checkRowIndex(this, row);
			}
			for (final int column : selectedColumns) {
				MatrixUtils.checkColumnIndex(this, column);
			}
		}
		return new BigMatrixImpl(subMatrixData, false);
	}

	/**
	 * Replace the submatrix starting at <code>row, column</code> using data in
	 * the input <code>subMatrix</code> array. Indexes are 0-based.
	 * <p>
	 * Example:<br>
	 * Starting with
	 * 
	 * <pre>
	 * 1  2  3  4
	 * 5  6  7  8
	 * 9  0  1  2
	 * </pre>
	 * 
	 * and <code>subMatrix = {{3, 4} {5,6}}</code>, invoking
	 * <code>setSubMatrix(subMatrix,1,1))</code> will result in
	 * 
	 * <pre>
	 * 1  2  3  4
	 * 5  3  4  8
	 * 9  5  6  2
	 * </pre>
	 * </p>
	 * 
	 * @param subMatrix
	 *            array containing the submatrix replacement data
	 * @param row
	 *            row coordinate of the top, left element to be replaced
	 * @param column
	 *            column coordinate of the top, left element to be replaced
	 * @throws MatrixIndexException
	 *             if subMatrix does not fit into this matrix from element in
	 *             (row, column)
	 * @throws IllegalArgumentException
	 *             if <code>subMatrix</code> is not rectangular (not all rows
	 *             have the same length) or empty
	 * @throws NullPointerException
	 *             if <code>subMatrix</code> is null
	 * @since 1.1
	 */
	public void setSubMatrix(BigDecimal[][] subMatrix, int row, int column) throws MatrixIndexException {

		final int nRows = subMatrix.length;
		if (nRows == 0) {
			throw MathRuntimeException.createIllegalArgumentException("matrix must have at least one row");
		}

		final int nCols = subMatrix[0].length;
		if (nCols == 0) {
			throw MathRuntimeException.createIllegalArgumentException("matrix must have at least one column");
		}

		for (int r = 1; r < nRows; r++) {
			if (subMatrix[r].length != nCols) {
				throw MathRuntimeException.createIllegalArgumentException(
						"some rows have length {0} while others have length {1}", nCols, subMatrix[r].length);
			}
		}

		if (data == null) {
			if (row > 0) {
				throw MathRuntimeException.createIllegalStateException("first {0} rows are not initialized yet", row);
			}
			if (column > 0) {
				throw MathRuntimeException.createIllegalStateException("first {0} columns are not initialized yet",
						column);
			}
			data = new BigDecimal[nRows][nCols];
			System.arraycopy(subMatrix, 0, data, 0, subMatrix.length);
		} else {
			MatrixUtils.checkRowIndex(this, row);
			MatrixUtils.checkColumnIndex(this, column);
			MatrixUtils.checkRowIndex(this, nRows + row - 1);
			MatrixUtils.checkColumnIndex(this, nCols + column - 1);
		}
		for (int i = 0; i < nRows; i++) {
			System.arraycopy(subMatrix[i], 0, data[row + i], column, nCols);
		}

		lu = null;

	}

	/**
	 * Returns the entries in row number <code>row</code> as a row matrix. Row
	 * indices start at 0.
	 *
	 * @param row
	 *            the row to be fetched
	 * @return row matrix
	 * @throws MatrixIndexException
	 *             if the specified row index is invalid
	 */
	public BigMatrix getRowMatrix(int row) throws MatrixIndexException {
		MatrixUtils.checkRowIndex(this, row);
		final int ncols = this.getColumnDimension();
		final BigDecimal[][] out = new BigDecimal[1][ncols];
		System.arraycopy(data[row], 0, out[0], 0, ncols);
		return new BigMatrixImpl(out, false);
	}

	/**
	 * Returns the entries in column number <code>column</code> as a column
	 * matrix. Column indices start at 0.
	 *
	 * @param column
	 *            the column to be fetched
	 * @return column matrix
	 * @throws MatrixIndexException
	 *             if the specified column index is invalid
	 */
	public BigMatrix getColumnMatrix(int column) throws MatrixIndexException {
		MatrixUtils.checkColumnIndex(this, column);
		final int nRows = this.getRowDimension();
		final BigDecimal[][] out = new BigDecimal[nRows][1];
		for (int row = 0; row < nRows; row++) {
			out[row][0] = data[row][column];
		}
		return new BigMatrixImpl(out, false);
	}

	/**
	 * Returns the entries in row number <code>row</code> as an array.
	 * <p>
	 * Row indices start at 0. A <code>MatrixIndexException</code> is thrown
	 * unless <code>0 <= row < rowDimension.</code>
	 * </p>
	 *
	 * @param row
	 *            the row to be fetched
	 * @return array of entries in the row
	 * @throws MatrixIndexException
	 *             if the specified row index is not valid
	 */
	public BigDecimal[] getRow(int row) throws MatrixIndexException {
		MatrixUtils.checkRowIndex(this, row);
		final int ncols = this.getColumnDimension();
		final BigDecimal[] out = new BigDecimal[ncols];
		System.arraycopy(data[row], 0, out, 0, ncols);
		return out;
	}

	/**
	 * Returns the entries in row number <code>row</code> as an array of double
	 * values.
	 * <p>
	 * Row indices start at 0. A <code>MatrixIndexException</code> is thrown
	 * unless <code>0 <= row < rowDimension.</code>
	 * </p>
	 *
	 * @param row
	 *            the row to be fetched
	 * @return array of entries in the row
	 * @throws MatrixIndexException
	 *             if the specified row index is not valid
	 */
	public double[] getRowAsDoubleArray(int row) throws MatrixIndexException {
		MatrixUtils.checkRowIndex(this, row);
		final int ncols = this.getColumnDimension();
		final double[] out = new double[ncols];
		for (int i = 0; i < ncols; i++) {
			out[i] = data[row][i].doubleValue();
		}
		return out;
	}

	/**
	 * Returns the entries in column number <code>col</code> as an array.
	 * <p>
	 * Column indices start at 0. A <code>MatrixIndexException</code> is thrown
	 * unless <code>0 <= column < columnDimension.</code>
	 * </p>
	 *
	 * @param col
	 *            the column to be fetched
	 * @return array of entries in the column
	 * @throws MatrixIndexException
	 *             if the specified column index is not valid
	 */
	public BigDecimal[] getColumn(int col) throws MatrixIndexException {
		MatrixUtils.checkColumnIndex(this, col);
		final int nRows = this.getRowDimension();
		final BigDecimal[] out = new BigDecimal[nRows];
		for (int i = 0; i < nRows; i++) {
			out[i] = data[i][col];
		}
		return out;
	}

	/**
	 * Returns the entries in column number <code>col</code> as an array of
	 * double values.
	 * <p>
	 * Column indices start at 0. A <code>MatrixIndexException</code> is thrown
	 * unless <code>0 <= column < columnDimension.</code>
	 * </p>
	 *
	 * @param col
	 *            the column to be fetched
	 * @return array of entries in the column
	 * @throws MatrixIndexException
	 *             if the specified column index is not valid
	 */
	public double[] getColumnAsDoubleArray(int col) throws MatrixIndexException {
		MatrixUtils.checkColumnIndex(this, col);
		final int nrows = this.getRowDimension();
		final double[] out = new double[nrows];
		for (int i = 0; i < nrows; i++) {
			out[i] = data[i][col].doubleValue();
		}
		return out;
	}

	/**
	 * Returns the entry in the specified row and column.
	 * <p>
	 * Row and column indices start at 0 and must satisfy
	 * <ul>
	 * <li><code>0 <= row < rowDimension</code></li>
	 * <li><code> 0 <= column < columnDimension</code></li>
	 * </ul>
	 * otherwise a <code>MatrixIndexException</code> is thrown.
	 * </p>
	 *
	 * @param row
	 *            row location of entry to be fetched
	 * @param column
	 *            column location of entry to be fetched
	 * @return matrix entry in row,column
	 * @throws MatrixIndexException
	 *             if the row or column index is not valid
	 */
	public BigDecimal getEntry(int row, int column) throws MatrixIndexException {
		try {
			return data[row][column];
		} catch (ArrayIndexOutOfBoundsException e) {
			throw new MatrixIndexException("no entry at indices ({0}, {1}) in a {2}x{3} matrix", row, column,
					getRowDimension(), getColumnDimension());
		}
	}

	/**
	 * Returns the entry in the specified row and column as a double.
	 * <p>
	 * Row and column indices start at 0 and must satisfy
	 * <ul>
	 * <li><code>0 <= row < rowDimension</code></li>
	 * <li><code> 0 <= column < columnDimension</code></li>
	 * </ul>
	 * otherwise a <code>MatrixIndexException</code> is thrown.
	 * </p>
	 *
	 * @param row
	 *            row location of entry to be fetched
	 * @param column
	 *            column location of entry to be fetched
	 * @return matrix entry in row,column
	 * @throws MatrixIndexException
	 *             if the row or column index is not valid
	 */
	public double getEntryAsDouble(int row, int column) throws MatrixIndexException {
		return getEntry(row, column).doubleValue();
	}

	/**
	 * Returns the transpose matrix.
	 *
	 * @return transpose matrix
	 */
	public BigMatrix transpose() {
		final int nRows = this.getRowDimension();
		final int nCols = this.getColumnDimension();
		final BigDecimal[][] outData = new BigDecimal[nCols][nRows];
		for (int row = 0; row < nRows; row++) {
			final BigDecimal[] dataRow = data[row];
			for (int col = 0; col < nCols; col++) {
				outData[col][row] = dataRow[col];
			}
		}
		return new BigMatrixImpl(outData, false);
	}

	/**
	 * Returns the inverse matrix if this matrix is invertible.
	 * 
	 * @return inverse matrix
	 * @throws InvalidMatrixException
	 *             if this is not invertible
	 */
	public BigMatrix inverse() throws InvalidMatrixException {
		return solve(MatrixUtils.createBigIdentityMatrix(getRowDimension()));
	}

	/**
	 * Returns the determinant of this matrix.
	 *
	 * @return determinant
	 * @throws InvalidMatrixException
	 *             if matrix is not square
	 */
	public BigDecimal getDeterminant() throws InvalidMatrixException {
		if (!isSquare()) {
			throw new NonSquareMatrixException(getRowDimension(), getColumnDimension());
		}
		if (isSingular()) { // note: this has side effect of attempting LU
							// decomp if lu == null
			return ZERO;
		} else {
			BigDecimal det = (parity == 1) ? ONE : ONE.negate();
			for (int i = 0; i < this.getRowDimension(); i++) {
				det = det.multiply(lu[i][i]);
			}
			return det;
		}
	}

	/**
	 * Is this a square matrix?
	 * 
	 * @return true if the matrix is square (rowDimension = columnDimension)
	 */
	public boolean isSquare() {
		return (this.getColumnDimension() == this.getRowDimension());
	}

	/**
	 * Is this a singular matrix?
	 * 
	 * @return true if the matrix is singular
	 */
	public boolean isSingular() {
		if (lu == null) {
			try {
				luDecompose();
				return false;
			} catch (InvalidMatrixException ex) {
				return true;
			}
		} else { // LU decomp must have been successfully performed
			return false; // so the matrix is not singular
		}
	}

	/**
	 * Returns the number of rows in the matrix.
	 *
	 * @return rowDimension
	 */
	public int getRowDimension() {
		return data.length;
	}

	/**
	 * Returns the number of columns in the matrix.
	 *
	 * @return columnDimension
	 */
	public int getColumnDimension() {
		return data[0].length;
	}

	/**
	 * Returns the <a href="http://mathworld.wolfram.com/MatrixTrace.html">
	 * trace</a> of the matrix (the sum of the elements on the main diagonal).
	 *
	 * @return trace
	 * 
	 * @throws IllegalArgumentException
	 *             if this matrix is not square.
	 */
	public BigDecimal getTrace() throws IllegalArgumentException {
		if (!isSquare()) {
			throw new NonSquareMatrixException(getRowDimension(), getColumnDimension());
		}
		BigDecimal trace = data[0][0];
		for (int i = 1; i < this.getRowDimension(); i++) {
			trace = trace.add(data[i][i]);
		}
		return trace;
	}

	/**
	 * Returns the result of multiplying this by the vector <code>v</code>.
	 *
	 * @param v
	 *            the vector to operate on
	 * @return this*v
	 * @throws IllegalArgumentException
	 *             if columnDimension != v.size()
	 */
	public BigDecimal[] operate(BigDecimal[] v) throws IllegalArgumentException {
		if (v.length != getColumnDimension()) {
			throw MathRuntimeException.createIllegalArgumentException(
					"vector length mismatch: got {0} but expected {1}", v.length, getColumnDimension());
		}
		final int nRows = this.getRowDimension();
		final int nCols = this.getColumnDimension();
		final BigDecimal[] out = new BigDecimal[nRows];
		for (int row = 0; row < nRows; row++) {
			BigDecimal sum = ZERO;
			for (int i = 0; i < nCols; i++) {
				sum = sum.add(data[row][i].multiply(v[i]));
			}
			out[row] = sum;
		}
		return out;
	}

	/**
	 * Returns the result of multiplying this by the vector <code>v</code>.
	 *
	 * @param v
	 *            the vector to operate on
	 * @return this*v
	 * @throws IllegalArgumentException
	 *             if columnDimension != v.size()
	 */
	public BigDecimal[] operate(double[] v) throws IllegalArgumentException {
		final BigDecimal bd[] = new BigDecimal[v.length];
		for (int i = 0; i < bd.length; i++) {
			bd[i] = new BigDecimal(v[i]);
		}
		return operate(bd);
	}

	/**
	 * Returns the (row) vector result of premultiplying this by the vector
	 * <code>v</code>.
	 *
	 * @param v
	 *            the row vector to premultiply by
	 * @return v*this
	 * @throws IllegalArgumentException
	 *             if rowDimension != v.size()
	 */
	public BigDecimal[] preMultiply(BigDecimal[] v) throws IllegalArgumentException {
		final int nRows = this.getRowDimension();
		if (v.length != nRows) {
			throw MathRuntimeException.createIllegalArgumentException(
					"vector length mismatch: got {0} but expected {1}", v.length, nRows);
		}
		final int nCols = this.getColumnDimension();
		final BigDecimal[] out = new BigDecimal[nCols];
		for (int col = 0; col < nCols; col++) {
			BigDecimal sum = ZERO;
			for (int i = 0; i < nRows; i++) {
				sum = sum.add(data[i][col].multiply(v[i]));
			}
			out[col] = sum;
		}
		return out;
	}

	/**
	 * Returns a matrix of (column) solution vectors for linear systems with
	 * coefficient matrix = this and constant vectors = columns of
	 * <code>b</code>.
	 *
	 * @param b
	 *            array of constants forming RHS of linear systems to to solve
	 * @return solution array
	 * @throws IllegalArgumentException
	 *             if this.rowDimension != row dimension
	 * @throws InvalidMatrixException
	 *             if this matrix is not square or is singular
	 */
	public BigDecimal[] solve(BigDecimal[] b) throws IllegalArgumentException, InvalidMatrixException {
		final int nRows = this.getRowDimension();
		if (b.length != nRows) {
			throw MathRuntimeException.createIllegalArgumentException(
					"vector length mismatch: got {0} but expected {1}", b.length, nRows);
		}
		final BigMatrix bMatrix = new BigMatrixImpl(b);
		final BigDecimal[][] solution = ((BigMatrixImpl) (solve(bMatrix))).getDataRef();
		final BigDecimal[] out = new BigDecimal[nRows];
		for (int row = 0; row < nRows; row++) {
			out[row] = solution[row][0];
		}
		return out;
	}

	/**
	 * Returns a matrix of (column) solution vectors for linear systems with
	 * coefficient matrix = this and constant vectors = columns of
	 * <code>b</code>.
	 *
	 * @param b
	 *            array of constants forming RHS of linear systems to to solve
	 * @return solution array
	 * @throws IllegalArgumentException
	 *             if this.rowDimension != row dimension
	 * @throws InvalidMatrixException
	 *             if this matrix is not square or is singular
	 */
	public BigDecimal[] solve(double[] b) throws IllegalArgumentException, InvalidMatrixException {
		final BigDecimal bd[] = new BigDecimal[b.length];
		for (int i = 0; i < bd.length; i++) {
			bd[i] = new BigDecimal(b[i]);
		}
		return solve(bd);
	}

	/**
	 * Returns a matrix of (column) solution vectors for linear systems with
	 * coefficient matrix = this and constant vectors = columns of
	 * <code>b</code>.
	 *
	 * @param b
	 *            matrix of constant vectors forming RHS of linear systems to to
	 *            solve
	 * @return matrix of solution vectors
	 * @throws IllegalArgumentException
	 *             if this.rowDimension != row dimension
	 * @throws InvalidMatrixException
	 *             if this matrix is not square or is singular
	 */
	public BigMatrix solve(BigMatrix b) throws IllegalArgumentException, InvalidMatrixException {
		if (b.getRowDimension() != getRowDimension()) {
			throw MathRuntimeException.createIllegalArgumentException(
					"dimensions mismatch: got {0}x{1} but expected {2}x{3}", b.getRowDimension(),
					b.getColumnDimension(), getRowDimension(), "n");
		}
		if (!isSquare()) {
			throw new NonSquareMatrixException(getRowDimension(), getColumnDimension());
		}
		if (this.isSingular()) { // side effect: compute LU decomp
			throw new SingularMatrixException();
		}

		final int nCol = this.getColumnDimension();
		final int nColB = b.getColumnDimension();
		final int nRowB = b.getRowDimension();

		// Apply permutations to b
		final BigDecimal[][] bp = new BigDecimal[nRowB][nColB];
		for (int row = 0; row < nRowB; row++) {
			final BigDecimal[] bpRow = bp[row];
			for (int col = 0; col < nColB; col++) {
				bpRow[col] = b.getEntry(permutation[row], col);
			}
		}

		// Solve LY = b
		for (int col = 0; col < nCol; col++) {
			for (int i = col + 1; i < nCol; i++) {
				final BigDecimal[] bpI = bp[i];
				final BigDecimal[] luI = lu[i];
				for (int j = 0; j < nColB; j++) {
					bpI[j] = bpI[j].subtract(bp[col][j].multiply(luI[col]));
				}
			}
		}

		// Solve UX = Y
		for (int col = nCol - 1; col >= 0; col--) {
			final BigDecimal[] bpCol = bp[col];
			final BigDecimal luDiag = lu[col][col];
			for (int j = 0; j < nColB; j++) {
				bpCol[j] = bpCol[j].divide(luDiag, scale, roundingMode);
			}
			for (int i = 0; i < col; i++) {
				final BigDecimal[] bpI = bp[i];
				final BigDecimal[] luI = lu[i];
				for (int j = 0; j < nColB; j++) {
					bpI[j] = bpI[j].subtract(bp[col][j].multiply(luI[col]));
				}
			}
		}

		return new BigMatrixImpl(bp, false);

	}

	/**
	 * Computes a new <a href=
	 * "http://www.math.gatech.edu/~bourbaki/math2601/Web-notes/2num.pdf"> LU
	 * decompostion</a> for this matrix, storing the result for use by other
	 * methods.
	 * <p>
	 * <strong>Implementation Note</strong>:<br>
	 * Uses <a href=
	 * "http://www.damtp.cam.ac.uk/user/fdl/people/sd/lectures/nummeth98/linear.htm">
	 * Crout's algortithm</a>, with partial pivoting.
	 * </p>
	 * <p>
	 * <strong>Usage Note</strong>:<br>
	 * This method should rarely be invoked directly. Its only use is to force
	 * recomputation of the LU decomposition when changes have been made to the
	 * underlying data using direct array references. Changes made using setXxx
	 * methods will trigger recomputation when needed automatically.
	 * </p>
	 *
	 * @throws InvalidMatrixException
	 *             if the matrix is non-square or singular.
	 */
	public void luDecompose() throws InvalidMatrixException {

		final int nRows = this.getRowDimension();
		final int nCols = this.getColumnDimension();
		if (nRows != nCols) {
			throw new NonSquareMatrixException(getRowDimension(), getColumnDimension());
		}
		lu = this.getData();

		// Initialize permutation array and parity
		permutation = new int[nRows];
		for (int row = 0; row < nRows; row++) {
			permutation[row] = row;
		}
		parity = 1;

		// Loop over columns
		for (int col = 0; col < nCols; col++) {

			BigDecimal sum = ZERO;

			// upper
			for (int row = 0; row < col; row++) {
				final BigDecimal[] luRow = lu[row];
				sum = luRow[col];
				for (int i = 0; i < row; i++) {
					sum = sum.subtract(luRow[i].multiply(lu[i][col]));
				}
				luRow[col] = sum;
			}

			// lower
			int max = col; // permutation row
			BigDecimal largest = ZERO;
			for (int row = col; row < nRows; row++) {
				final BigDecimal[] luRow = lu[row];
				sum = luRow[col];
				for (int i = 0; i < col; i++) {
					sum = sum.subtract(luRow[i].multiply(lu[i][col]));
				}
				luRow[col] = sum;

				// maintain best permutation choice
				if (sum.abs().compareTo(largest) == 1) {
					largest = sum.abs();
					max = row;
				}
			}

			// Singularity check
			if (lu[max][col].abs().compareTo(TOO_SMALL) <= 0) {
				lu = null;
				throw new SingularMatrixException();
			}

			// Pivot if necessary
			if (max != col) {
				BigDecimal tmp = ZERO;
				for (int i = 0; i < nCols; i++) {
					tmp = lu[max][i];
					lu[max][i] = lu[col][i];
					lu[col][i] = tmp;
				}
				int temp = permutation[max];
				permutation[max] = permutation[col];
				permutation[col] = temp;
				parity = -parity;
			}

			// Divide the lower elements by the "winning" diagonal elt.
			final BigDecimal luDiag = lu[col][col];
			for (int row = col + 1; row < nRows; row++) {
				final BigDecimal[] luRow = lu[row];
				luRow[col] = luRow[col].divide(luDiag, scale, roundingMode);
			}

		}

	}

	/**
	 * Get a string representation for this matrix.
	 * 
	 * @return a string representation for this matrix
	 */
	@Override
	public String toString() {
		StringBuffer res = new StringBuffer();
		res.append("BigMatrixImpl{");
		if (data != null) {
			for (int i = 0; i < data.length; i++) {
				if (i > 0) {
					res.append(",");
				}
				res.append("{");
				for (int j = 0; j < data[0].length; j++) {
					if (j > 0) {
						res.append(",");
					}
					res.append(data[i][j]);
				}
				res.append("}");
			}
		}
		res.append("}");
		return res.toString();
	}

	/**
	 * Returns true iff <code>object</code> is a <code>BigMatrixImpl</code>
	 * instance with the same dimensions as this and all corresponding matrix
	 * entries are equal. BigDecimal.equals is used to compare corresponding
	 * entries.
	 * 
	 * @param object
	 *            the object to test equality against.
	 * @return true if object equals this
	 */
	@Override
	public boolean equals(Object object) {
		if (object == this) {
			return true;
		}
		if (object instanceof BigMatrixImpl == false) {
			return false;
		}
		final BigMatrix m = (BigMatrix) object;
		final int nRows = getRowDimension();
		final int nCols = getColumnDimension();
		if (m.getColumnDimension() != nCols || m.getRowDimension() != nRows) {
			return false;
		}
		for (int row = 0; row < nRows; row++) {
			final BigDecimal[] dataRow = data[row];
			for (int col = 0; col < nCols; col++) {
				if (!dataRow[col].equals(m.getEntry(row, col))) {
					return false;
				}
			}
		}
		return true;
	}

	/**
	 * Computes a hashcode for the matrix.
	 * 
	 * @return hashcode for matrix
	 */
	@Override
	public int hashCode() {
		int ret = 7;
		final int nRows = getRowDimension();
		final int nCols = getColumnDimension();
		ret = ret * 31 + nRows;
		ret = ret * 31 + nCols;
		for (int row = 0; row < nRows; row++) {
			final BigDecimal[] dataRow = data[row];
			for (int col = 0; col < nCols; col++) {
				ret = ret * 31 + (11 * (row + 1) + 17 * (col + 1)) * dataRow[col].hashCode();
			}
		}
		return ret;
	}

	// ------------------------ Protected methods

	/**
	 * Returns the LU decomposition as a BigMatrix. Returns a fresh copy of the
	 * cached LU matrix if this has been computed; otherwise the composition is
	 * computed and cached for use by other methods. Since a copy is returned in
	 * either case, changes to the returned matrix do not affect the LU
	 * decomposition property.
	 * <p>
	 * The matrix returned is a compact representation of the LU decomposition.
	 * Elements below the main diagonal correspond to entries of the "L" matrix;
	 * elements on and above the main diagonal correspond to entries of the "U"
	 * matrix.
	 * </p>
	 * <p>
	 * Example:
	 * 
	 * <pre>
	 * 
	 *     Returned matrix                L                  U
	 *         2  3  1                   1  0  0            2  3  1          
	 *         5  4  6                   5  1  0            0  4  6
	 *         1  7  8                   1  7  1            0  0  8
	 * </pre>
	 * 
	 * The L and U matrices satisfy the matrix equation LU = permuteRows(this),
	 * <br>
	 * where permuteRows reorders the rows of the matrix to follow the order
	 * determined by the <a href=#getPermutation()>permutation</a> property.
	 * </p>
	 * 
	 * @return LU decomposition matrix
	 * @throws InvalidMatrixException
	 *             if the matrix is non-square or singular.
	 */
	protected BigMatrix getLUMatrix() throws InvalidMatrixException {
		if (lu == null) {
			luDecompose();
		}
		return new BigMatrixImpl(lu);
	}

	/**
	 * Returns the permutation associated with the lu decomposition. The entries
	 * of the array represent a permutation of the numbers 0, ... , nRows - 1.
	 * <p>
	 * Example: permutation = [1, 2, 0] means current 2nd row is first, current
	 * third row is second and current first row is last.
	 * </p>
	 * <p>
	 * Returns a fresh copy of the array.
	 * </p>
	 * 
	 * @return the permutation
	 */
	protected int[] getPermutation() {
		final int[] out = new int[permutation.length];
		System.arraycopy(permutation, 0, out, 0, permutation.length);
		return out;
	}

	// ------------------------ Private methods

	/**
	 * Returns a fresh copy of the underlying data array.
	 *
	 * @return a copy of the underlying data array.
	 */
	private BigDecimal[][] copyOut() {
		final int nRows = this.getRowDimension();
		final BigDecimal[][] out = new BigDecimal[nRows][this.getColumnDimension()];
		// can't copy 2-d array in one shot, otherwise get row references
		for (int i = 0; i < nRows; i++) {
			System.arraycopy(data[i], 0, out[i], 0, data[i].length);
		}
		return out;
	}

	/**
	 * Replaces data with a fresh copy of the input array.
	 * <p>
	 * Verifies that the input array is rectangular and non-empty.
	 * </p>
	 *
	 * @param in
	 *            data to copy in
	 * @throws IllegalArgumentException
	 *             if input array is emtpy or not rectangular
	 * @throws NullPointerException
	 *             if input array is null
	 */
	private void copyIn(BigDecimal[][] in) {
		setSubMatrix(in, 0, 0);
	}

	/**
	 * Replaces data with a fresh copy of the input array.
	 *
	 * @param in
	 *            data to copy in
	 */
	private void copyIn(double[][] in) {
		final int nRows = in.length;
		final int nCols = in[0].length;
		data = new BigDecimal[nRows][nCols];
		for (int i = 0; i < nRows; i++) {
			final BigDecimal[] dataI = data[i];
			final double[] inI = in[i];
			for (int j = 0; j < nCols; j++) {
				dataI[j] = new BigDecimal(inI[j]);
			}
		}
		lu = null;
	}

	/**
	 * Replaces data with BigDecimals represented by the strings in the input
	 * array.
	 *
	 * @param in
	 *            data to copy in
	 */
	private void copyIn(String[][] in) {
		final int nRows = in.length;
		final int nCols = in[0].length;
		data = new BigDecimal[nRows][nCols];
		for (int i = 0; i < nRows; i++) {
			final BigDecimal[] dataI = data[i];
			final String[] inI = in[i];
			for (int j = 0; j < nCols; j++) {
				dataI[j] = new BigDecimal(inI[j]);
			}
		}
		lu = null;
	}

}
