densematrix.java

另一个功能更强大的矩阵运算软件开源代码

/*
 * Copyright (C) 2003-2006 Bjørn-Ove Heimsund
 * 
 * This file is part of MTJ.
 * 
 * This library is free software; you can redistribute it and/or modify it
 * under the terms of the GNU Lesser General Public License as published by the
 * Free Software Foundation; either version 2.1 of the License, or (at your
 * option) any later version.
 * 
 * This library is distributed in the hope that it will be useful, but WITHOUT
 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License
 * for more details.
 * 
 * You should have received a copy of the GNU Lesser General Public License
 * along with this library; if not, write to the Free Software Foundation,
 * Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
 */

package no.uib.cipr.matrix;

import java.io.IOException;

import no.uib.cipr.matrix.io.MatrixInfo;
import no.uib.cipr.matrix.io.MatrixSize;
import no.uib.cipr.matrix.io.MatrixVectorReader;

import org.netlib.blas.BLAS;
import org.netlib.lapack.LAPACK;
import org.netlib.util.intW;

/**
 * Dense matrix. It is a good all-round matrix structure, with fast access and
 * efficient algebraic operations. The matrix
 * <p>
 * <table border="1">
 * <tr>
 * <td>a<sub>11</sub></td>
 * <td>a<sub>12</sub></td>
 * <td>a<sub>13</sub></td>
 * <td>a<sub>14</sub></td>
 * </tr>
 * <tr>
 * <td>a<sub>21</sub></td>
 * <td>a<sub>22</sub></td>
 * <td>a<sub>23</sub></td>
 * <td>a<sub>24</sub></td>
 * </tr>
 * <tr>
 * <td>a<sub>31</sub></td>
 * <td>a<sub>32</sub></td>
 * <td>a<sub>33</sub></td>
 * <td>a<sub>34</sub></td>
 * </tr>
 * <tr>
 * <td>a<sub>41</sub></td>
 * <td>a<sub>42</sub></td>
 * <td>a<sub>43</sub></td>
 * <td>a<sub>44</sub></td>
 * </tr>
 * </table>
 * </p>
 * <p>
 * is stored column major in a single array, as follows:
 * </p>
 * <p>
 * <table border="1">
 * <tr>
 * <td>a<sub>11</sub></td>
 * <td>a<sub>21</sub></td>
 * <td>a<sub>31</sub></td>
 * <td>a<sub>41</sub></td>
 * <td>a<sub>12</sub></td>
 * <td>a<sub>22</sub></td>
 * <td>a<sub>32</sub></td>
 * <td>a<sub>42</sub></td>
 * <td>a<sub>13</sub></td>
 * <td>a<sub>23</sub></td>
 * <td>a<sub>33</sub></td>
 * <td>a<sub>43</sub></td>
 * <td>a<sub>14</sub></td>
 * <td>a<sub>24</sub></td>
 * <td>a<sub>34</sub></td>
 * <td>a<sub>44</sub></td>
 * </tr>
 * </table>
 * </p>
 */
public class DenseMatrix extends AbstractDenseMatrix {

    /**
     * Constructor for DenseMatrix
     * 
     * @param r
     *            Reader to get the matrix from
     */
    public DenseMatrix(MatrixVectorReader r) throws IOException {
        // Start with a zero-sized matrix
        super(0, 0);

        // Get matrix information. Use the header if present, else use a safe
        // default
        MatrixInfo info = null;
        if (r.hasInfo())
            info = r.readMatrixInfo();
        else
            info = new MatrixInfo(true, MatrixInfo.MatrixField.Real,
                    MatrixInfo.MatrixSymmetry.General);
        MatrixSize size = r.readMatrixSize(info);

        // Resize the matrix to correct size
        numRows = size.numRows();
        numColumns = size.numColumns();
        data = new double[numRows * numColumns];

        // Check that the matrix is in an acceptable format
        if (info.isPattern())
            throw new UnsupportedOperationException(
                    "Pattern matrices are not supported");
        if (info.isComplex())
            throw new UnsupportedOperationException(
                    "Complex matrices are not supported");

        // Read the entries, in either coordinate or array format
        if (info.isCoordinate()) {

            // Read coordinate data
            int nz = size.numEntries();
            int[] row = new int[nz];
            int[] column = new int[nz];
            double[] entry = new double[nz];
            r.readCoordinate(row, column, entry);

            // Shift indices from 1-offset to 0-offset
            r.add(-1, row);
            r.add(-1, column);

            // Store them
            for (int i = 0; i < nz; ++i)
                set(row[i], column[i], entry[i]);

        } else
            // info.isArray()
            r.readArray(data);

        // Put in missing entries from symmetry or skew symmetry
        if (info.isSymmetric())
            for (int i = 0; i < numRows; ++i)
                for (int j = 0; j < i; ++j)
                    set(j, i, get(i, j));
        else if (info.isSkewSymmetric())
            for (int i = 0; i < numRows; ++i)
                for (int j = 0; j < i; ++j)
                    set(j, i, -get(i, j));
    }

    /**
     * Constructor for DenseMatrix
     * 
     * @param numRows
     *            Number of rows
     * @param numColumns
     *            Number of columns
     */
    public DenseMatrix(int numRows, int numColumns) {
        super(numRows, numColumns);
    }

    /**
     * Constructor for DenseMatrix
     * 
     * @param A
     *            Matrix to copy. A deep copy is made
     */
    public DenseMatrix(Matrix A) {
        super(A);
    }

    /**
     * Constructor for DenseMatrix
     * 
     * @param A
     *            Matrix to copy contents from
     * @param deep
     *            If true, <code>A</code> is copied, else a shallow copy is
     *            made and the matrices share underlying storage. For this,
     *            <code>A</code> must be a dense matrix
     */
    public DenseMatrix(Matrix A, boolean deep) {
        super(A, deep);
    }

    /**
     * Constructor for DenseMatrix. Builds the matrix from a vector
     * 
     * @param x
     *            Vector to copy from. This will form this matrix' single column
     * @param deep
     *            If true, x is copied, if false, the internal storage of this
     *            matrix is the same as that of the vector. In that case,
     *            <code>x</code> must be a <code>DenseVector</code>
     */
    public DenseMatrix(Vector x, boolean deep) {
        super(x.size(), 1);

        if (deep)
            for (VectorEntry e : x)
                set(e.index(), 0, e.get());
        else {
            if (!(x instanceof DenseVector))
                throw new IllegalArgumentException("x must be a DenseVector");
            data = ((DenseVector) x).getData();
        }
    }

    /**
     * Constructor for DenseMatrix. Builds the matrix from a vector
     * 
     * @param x
     *            The vector which forms this matrix' single column. It is
     *            copied, not referenced
     */
    public DenseMatrix(Vector x) {
        this(x, true);
    }

    /**
     * Constructor for DenseMatrix. Builds the matrix from vectors. Each vector
     * will correspond to a column of the matrix
     * 
     * @param x
     *            Vectors which forms the columns of this matrix. Every vector
     *            must have the same size
     */
    public DenseMatrix(Vector[] x) {
        super(x[0].size(), x.length);

        // Ensure correct sizes
        for (Vector v : x)
            if (v.size() != numRows)
                throw new IllegalArgumentException(
                        "All vectors must be of the same size");

        // Copy the contents
        for (int j = 0; j < x.length; ++j)
            for (VectorEntry e : x[j])
                set(e.index(), j, e.get());
    }

    /**
     * Constructor for DenseMatrix. Copies from the passed array
     * 
     * @param values
     *            Arrays to copy from. Every sub-array must have the same size
     */
    public DenseMatrix(double[][] values) {
        super(values.length, values[0].length);