001    /*
002     * Copyright (C) 2008-2010 by Holger Arndt
003     *
004     * This file is part of the Universal Java Matrix Package (UJMP).
005     * See the NOTICE file distributed with this work for additional
006     * information regarding copyright ownership and licensing.
007     *
008     * UJMP is free software; you can redistribute it and/or modify
009     * it under the terms of the GNU Lesser General Public License as
010     * published by the Free Software Foundation; either version 2
011     * of the License, or (at your option) any later version.
012     *
013     * UJMP is distributed in the hope that it will be useful,
014     * but WITHOUT ANY WARRANTY; without even the implied warranty of
015     * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
016     * GNU Lesser General Public License for more details.
017     *
018     * You should have received a copy of the GNU Lesser General Public
019     * License along with UJMP; if not, write to the
020     * Free Software Foundation, Inc., 51 Franklin St, Fifth Floor,
021     * Boston, MA  02110-1301  USA
022     */
023    
024    package org.ujmp.ejml;
025    
026    import java.io.IOException;
027    import java.io.ObjectInputStream;
028    import java.io.ObjectOutputStream;
029    
030    import org.ejml.alg.dense.decomposition.CholeskyDecomposition;
031    import org.ejml.alg.dense.decomposition.DecompositionFactory;
032    import org.ejml.alg.dense.decomposition.EigenDecomposition;
033    import org.ejml.alg.dense.decomposition.LUDecomposition;
034    import org.ejml.data.Complex64F;
035    import org.ejml.data.DenseMatrix64F;
036    import org.ejml.ops.CommonOps;
037    import org.ujmp.core.Matrix;
038    import org.ujmp.core.doublematrix.DenseDoubleMatrix2D;
039    import org.ujmp.core.doublematrix.stub.AbstractDenseDoubleMatrix2D;
040    import org.ujmp.core.exceptions.MatrixException;
041    import org.ujmp.core.interfaces.HasRowMajorDoubleArray1D;
042    import org.ujmp.core.interfaces.HasRowMajorDoubleArray2D;
043    import org.ujmp.core.interfaces.Wrapper;
044    import org.ujmp.ejml.calculation.Inv;
045    import org.ujmp.ejml.calculation.InvSPD;
046    import org.ujmp.ejml.calculation.QR;
047    import org.ujmp.ejml.calculation.SVD;
048    import org.ujmp.ejml.calculation.Solve;
049    
050    public class EJMLDenseDoubleMatrix2D extends AbstractDenseDoubleMatrix2D
051                    implements Wrapper<DenseMatrix64F>, HasRowMajorDoubleArray1D {
052            private static final long serialVersionUID = -3223474248020842822L;
053    
054            private transient DenseMatrix64F matrix;
055    
056            public EJMLDenseDoubleMatrix2D(long... size) {
057                    this.matrix = new DenseMatrix64F((int) size[ROW], (int) size[COLUMN]);
058            }
059    
060            public EJMLDenseDoubleMatrix2D(DenseMatrix64F m) {
061                    this.matrix = (DenseMatrix64F) m;
062            }
063    
064            public EJMLDenseDoubleMatrix2D(Matrix source) throws MatrixException {
065                    // if (source instanceof HasColumnMajorDoubleArray1D) {
066                    // final double[] data = ((HasColumnMajorDoubleArray1D) source)
067                    // .getColumnMajorDoubleArray1D();
068                    // this.matrix = new DenseMatrix64F((int) source.getRowCount(),
069                    // (int) source.getColumnCount(), false, data);
070                    // } else
071                    if (source instanceof HasRowMajorDoubleArray2D) {
072                            final double[][] data = ((HasRowMajorDoubleArray2D) source)
073                                            .getRowMajorDoubleArray2D();
074                            this.matrix = new DenseMatrix64F(data);
075                    } else if (source instanceof DenseDoubleMatrix2D) {
076                            this.matrix = new DenseMatrix64F((int) source.getRowCount(),
077                                            (int) source.getColumnCount());
078                            final DenseDoubleMatrix2D m2 = (DenseDoubleMatrix2D) source;
079                            for (int r = (int) source.getRowCount(); --r >= 0;) {
080                                    for (int c = (int) source.getColumnCount(); --c >= 0;) {
081                                            matrix.set(r, c, m2.getDouble(r, c));
082                                    }
083                            }
084                    } else {
085                            this.matrix = new DenseMatrix64F((int) source.getRowCount(),
086                                            (int) source.getColumnCount());
087                            for (long[] c : source.availableCoordinates()) {
088                                    setDouble(source.getAsDouble(c), c);
089                            }
090                    }
091            }
092    
093            public double getDouble(long row, long column) {
094                    return matrix.get((int) row, (int) column);
095            }
096    
097            public double getDouble(int row, int column) {
098                    return matrix.get(row, column);
099            }
100    
101            public long[] getSize() {
102                    return new long[] { matrix.numRows, matrix.numCols };
103            }
104    
105            public void setDouble(double value, long row, long column) {
106                    matrix.set((int) row, (int) column, value);
107            }
108    
109            public void setDouble(double value, int row, int column) {
110                    matrix.set(row, column, value);
111            }
112    
113            public DenseMatrix64F getWrappedObject() {
114                    return matrix;
115            }
116    
117            public void setWrappedObject(DenseMatrix64F object) {
118                    this.matrix = object;
119            }
120    
121            public Matrix transpose() {
122                    DenseMatrix64F ret = new DenseMatrix64F(matrix.numCols, matrix.numRows);
123                    CommonOps.transpose(matrix, ret);
124                    return new EJMLDenseDoubleMatrix2D(ret);
125            }
126    
127            public Matrix inv() {
128                    return Inv.INSTANCE.calc(this);
129            }
130    
131            public Matrix invSPD() {
132                    return InvSPD.INSTANCE.calc(this);
133            }
134    
135            public double det() {
136                    return CommonOps.det(matrix);
137            }
138    
139            public Matrix solve(Matrix b) {
140                    return Solve.INSTANCE.calc(this, b);
141            }
142    
143            public Matrix plus(double value) {
144                    DenseMatrix64F ret = new DenseMatrix64F(matrix.numRows, matrix.numCols);
145                    CommonOps.add(matrix, value, ret);
146                    return new EJMLDenseDoubleMatrix2D(ret);
147            }
148    
149            public Matrix plus(Matrix m) {
150                    if (m instanceof EJMLDenseDoubleMatrix2D) {
151                            DenseMatrix64F ret = new DenseMatrix64F(matrix.numRows,
152                                            matrix.numCols);
153                            CommonOps.add(matrix, ((EJMLDenseDoubleMatrix2D) m).matrix, ret);
154                            return new EJMLDenseDoubleMatrix2D(ret);
155                    } else {
156                            return super.plus(m);
157                    }
158            }
159    
160            public Matrix minus(Matrix m) {
161                    if (m instanceof EJMLDenseDoubleMatrix2D) {
162                            DenseMatrix64F ret = new DenseMatrix64F(matrix.numRows,
163                                            matrix.numCols);
164                            CommonOps.sub(matrix, ((EJMLDenseDoubleMatrix2D) m).matrix, ret);
165                            return new EJMLDenseDoubleMatrix2D(ret);
166                    } else {
167                            return super.minus(m);
168                    }
169            }
170    
171            public Matrix minus(double value) {
172                    DenseMatrix64F ret = new DenseMatrix64F(matrix.numRows, matrix.numCols);
173                    CommonOps.add(matrix, -value, ret);
174                    return new EJMLDenseDoubleMatrix2D(ret);
175            }
176    
177            public Matrix times(double value) {
178                    DenseMatrix64F ret = new DenseMatrix64F(matrix.numRows, matrix.numCols);
179                    CommonOps.scale(value, matrix, ret);
180                    return new EJMLDenseDoubleMatrix2D(ret);
181            }
182    
183            public Matrix divide(double value) {
184                    DenseMatrix64F ret = new DenseMatrix64F(matrix.numRows, matrix.numCols);
185                    CommonOps.scale(1.0 / value, matrix, ret);
186                    return new EJMLDenseDoubleMatrix2D(ret);
187            }
188    
189            public Matrix mtimes(Matrix m) {
190                    if (m instanceof EJMLDenseDoubleMatrix2D) {
191                            DenseMatrix64F b = ((EJMLDenseDoubleMatrix2D) m).getWrappedObject();
192                            DenseMatrix64F ret = new DenseMatrix64F(matrix.numRows, b.numCols);
193                            CommonOps.mult(matrix, b, ret);
194                            return new EJMLDenseDoubleMatrix2D(ret);
195                    } else {
196                            return super.mtimes(m);
197                    }
198            }
199    
200            public Matrix[] svd() {
201                    return SVD.INSTANCE.calc(this);
202            }
203    
204            public Matrix[] qr() {
205                    return QR.INSTANCE.calc(this);
206            }
207    
208            public Matrix chol() {
209                    CholeskyDecomposition chol = DecompositionFactory.chol();
210                    chol.decompose(matrix);
211                    Matrix l = new EJMLDenseDoubleMatrix2D(chol.getT(null));
212                    return l;
213            }
214    
215            public Matrix[] lu() {
216                    if (isSquare()) {
217                            LUDecomposition lu = DecompositionFactory.lu();
218                            lu.decompose(matrix);
219                            DenseMatrix64F lm = new DenseMatrix64F(matrix.numRows,
220                                            matrix.numCols);
221                            DenseMatrix64F um = new DenseMatrix64F(matrix.numRows,
222                                            matrix.numCols);
223                            lu.getLower(lm);
224                            lu.getUpper(um);
225    
226                            Matrix l = new EJMLDenseDoubleMatrix2D(lm);
227                            Matrix u = new EJMLDenseDoubleMatrix2D(um);
228                            Matrix p = new EJMLDenseDoubleMatrix2D(lu.getPivot(null));
229    
230                            return new Matrix[] { l, u, p };
231                    } else {
232                            return super.lu();
233                    }
234            }
235    
236            public Matrix[] eig() {
237                    EigenDecomposition eig = DecompositionFactory.eig();
238                    eig.decompose(matrix);
239    
240                    int N = matrix.numRows;
241    
242                    DenseMatrix64F D = new DenseMatrix64F(N, N);
243                    DenseMatrix64F V = new DenseMatrix64F(N, N);
244    
245                    for (int i = 0; i < N; i++) {
246                            Complex64F c = eig.getEigenvalue(i);
247    
248                            if (c.isReal()) {
249                                    D.set(i, i, c.real);
250    
251                                    DenseMatrix64F v = eig.getEigenVector(i);
252    
253                                    if (v != null) {
254                                            for (int j = 0; j < N; j++) {
255                                                    V.set(j, i, v.get(j, 0));
256                                            }
257                                    }
258                            }
259                    }
260    
261                    Matrix v = new EJMLDenseDoubleMatrix2D(V);
262                    Matrix d = new EJMLDenseDoubleMatrix2D(D);
263    
264                    return new Matrix[] { v, d };
265            }
266    
267            public Matrix copy() {
268                    Matrix m = new EJMLDenseDoubleMatrix2D(matrix.copy());
269                    if (getAnnotation() != null) {
270                            m.setAnnotation(getAnnotation().clone());
271                    }
272                    return m;
273            }
274    
275            private void readObject(ObjectInputStream s) throws IOException,
276                            ClassNotFoundException {
277                    s.defaultReadObject();
278                    int rows = (Integer) s.readObject();
279                    int columns = (Integer) s.readObject();
280                    double[] values = (double[]) s.readObject();
281                    matrix = new DenseMatrix64F(rows, columns, true, values);
282            }
283    
284            private void writeObject(ObjectOutputStream s) throws IOException,
285                            MatrixException {
286                    s.defaultWriteObject();
287                    s.writeObject(matrix.numRows);
288                    s.writeObject(matrix.numCols);
289                    s.writeObject(matrix.data);
290            }
291    
292            public double[] getRowMajorDoubleArray1D() {
293                    return matrix.getData();
294            }
295    
296    }