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.commonsmath; 025 026 import org.apache.commons.math.linear.CholeskyDecomposition; 027 import org.apache.commons.math.linear.CholeskyDecompositionImpl; 028 import org.apache.commons.math.linear.EigenDecomposition; 029 import org.apache.commons.math.linear.EigenDecompositionImpl; 030 import org.apache.commons.math.linear.LUDecomposition; 031 import org.apache.commons.math.linear.LUDecompositionImpl; 032 import org.apache.commons.math.linear.QRDecomposition; 033 import org.apache.commons.math.linear.QRDecompositionImpl; 034 import org.apache.commons.math.linear.RealMatrix; 035 import org.apache.commons.math.linear.SingularValueDecomposition; 036 import org.apache.commons.math.linear.SingularValueDecompositionImpl; 037 import org.apache.commons.math.util.MathUtils; 038 import org.ujmp.core.Coordinates; 039 import org.ujmp.core.Matrix; 040 import org.ujmp.core.doublematrix.stub.AbstractDenseDoubleMatrix2D; 041 import org.ujmp.core.exceptions.MatrixException; 042 import org.ujmp.core.interfaces.Wrapper; 043 044 public abstract class AbstractCommonsMathDenseDoubleMatrix2D extends 045 AbstractDenseDoubleMatrix2D implements Wrapper<RealMatrix> { 046 private static final long serialVersionUID = -1161807620507675926L; 047 048 private RealMatrix matrix = null; 049 050 public AbstractCommonsMathDenseDoubleMatrix2D(RealMatrix matrix) { 051 this.matrix = matrix; 052 } 053 054 public RealMatrix getWrappedObject() { 055 return matrix; 056 } 057 058 public void setWrappedObject(RealMatrix object) { 059 this.matrix = object; 060 } 061 062 public double getDouble(long row, long column) throws MatrixException { 063 return matrix.getEntry((int) row, (int) column); 064 } 065 066 public double getDouble(int row, int column) throws MatrixException { 067 return matrix.getEntry(row, column); 068 } 069 070 public void setDouble(double value, long row, long column) 071 throws MatrixException { 072 matrix.setEntry((int) row, (int) column, value); 073 } 074 075 public void setDouble(double value, int row, int column) 076 throws MatrixException { 077 matrix.setEntry((int) row, (int) column, value); 078 } 079 080 public long[] getSize() { 081 return matrix == null ? Coordinates.ZERO2D : new long[] { 082 matrix.getRowDimension(), matrix.getColumnDimension() }; 083 } 084 085 public Matrix transpose() { 086 return CommonsMathDenseDoubleMatrix2DFactory.INSTANCE.dense(matrix 087 .transpose()); 088 } 089 090 public Matrix inv() { 091 return CommonsMathDenseDoubleMatrix2DFactory.INSTANCE 092 .dense(new LUDecompositionImpl(matrix).getSolver().getInverse()); 093 } 094 095 public Matrix invSPD() { 096 try { 097 return CommonsMathDenseDoubleMatrix2DFactory.INSTANCE 098 .dense(new CholeskyDecompositionImpl(matrix).getSolver() 099 .getInverse()); 100 } catch (Exception e) { 101 throw new MatrixException(e); 102 } 103 } 104 105 public Matrix[] lu() { 106 LUDecomposition lu = new LUDecompositionImpl(matrix); 107 Matrix l = CommonsMathDenseDoubleMatrix2DFactory.INSTANCE.dense(lu 108 .getL()); 109 Matrix u = CommonsMathDenseDoubleMatrix2DFactory.INSTANCE.dense(lu 110 .getU()); 111 Matrix p = CommonsMathDenseDoubleMatrix2DFactory.INSTANCE.dense(lu 112 .getP()); 113 return new Matrix[] { l, u, p }; 114 } 115 116 public Matrix[] qr() { 117 QRDecomposition qr = new QRDecompositionImpl(matrix); 118 Matrix q = CommonsMathDenseDoubleMatrix2DFactory.INSTANCE.dense(qr 119 .getQ()); 120 Matrix r = CommonsMathDenseDoubleMatrix2DFactory.INSTANCE.dense(qr 121 .getR()); 122 return new Matrix[] { q, r }; 123 } 124 125 public Matrix[] svd() { 126 SingularValueDecomposition svd = new SingularValueDecompositionImpl( 127 matrix); 128 Matrix u = CommonsMathDenseDoubleMatrix2DFactory.INSTANCE.dense(svd 129 .getU()); 130 Matrix s = CommonsMathDenseDoubleMatrix2DFactory.INSTANCE.dense(svd 131 .getS()); 132 Matrix v = CommonsMathDenseDoubleMatrix2DFactory.INSTANCE.dense(svd 133 .getV()); 134 return new Matrix[] { u, s, v }; 135 } 136 137 public Matrix[] eig() { 138 EigenDecomposition evd = new EigenDecompositionImpl(matrix, 139 MathUtils.EPSILON); 140 Matrix v = CommonsMathDenseDoubleMatrix2DFactory.INSTANCE.dense(evd 141 .getV()); 142 Matrix d = CommonsMathDenseDoubleMatrix2DFactory.INSTANCE.dense(evd 143 .getD()); 144 return new Matrix[] { v, d }; 145 } 146 147 public Matrix chol() { 148 try { 149 CholeskyDecomposition chol = new CholeskyDecompositionImpl(matrix); 150 Matrix l = CommonsMathDenseDoubleMatrix2DFactory.INSTANCE 151 .dense(chol.getL()); 152 return l; 153 } catch (Exception e) { 154 throw new MatrixException(e); 155 } 156 } 157 158 public Matrix mtimes(Matrix m2) { 159 if (m2 instanceof AbstractCommonsMathDenseDoubleMatrix2D) { 160 return CommonsMathDenseDoubleMatrix2DFactory.INSTANCE 161 .dense(matrix 162 .multiply(((AbstractCommonsMathDenseDoubleMatrix2D) m2).matrix)); 163 } else { 164 return super.mtimes(m2); 165 } 166 } 167 168 public Matrix plus(Matrix m2) { 169 if (m2 instanceof AbstractCommonsMathDenseDoubleMatrix2D) { 170 return CommonsMathDenseDoubleMatrix2DFactory.INSTANCE.dense(matrix 171 .add(((AbstractCommonsMathDenseDoubleMatrix2D) m2).matrix)); 172 } else { 173 return super.plus(m2); 174 } 175 } 176 177 public Matrix minus(Matrix m2) { 178 if (m2 instanceof AbstractCommonsMathDenseDoubleMatrix2D) { 179 return CommonsMathDenseDoubleMatrix2DFactory.INSTANCE 180 .dense(matrix 181 .subtract(((AbstractCommonsMathDenseDoubleMatrix2D) m2).matrix)); 182 } else { 183 return super.minus(m2); 184 } 185 } 186 187 public Matrix times(double value) { 188 return CommonsMathDenseDoubleMatrix2DFactory.INSTANCE.dense(matrix 189 .scalarMultiply(value)); 190 } 191 192 public Matrix divide(double value) { 193 return CommonsMathDenseDoubleMatrix2DFactory.INSTANCE.dense(matrix 194 .scalarMultiply(1.0 / value)); 195 } 196 197 public Matrix plus(double value) { 198 return CommonsMathDenseDoubleMatrix2DFactory.INSTANCE.dense(matrix 199 .scalarAdd(value)); 200 } 201 202 public Matrix minus(double value) { 203 return CommonsMathDenseDoubleMatrix2DFactory.INSTANCE.dense(matrix 204 .scalarAdd(-value)); 205 } 206 207 public Matrix solve(Matrix b) { 208 if (b instanceof AbstractCommonsMathDenseDoubleMatrix2D) { 209 AbstractCommonsMathDenseDoubleMatrix2D b2 = (AbstractCommonsMathDenseDoubleMatrix2D) b; 210 if (isSquare()) { 211 RealMatrix ret = new LUDecompositionImpl(matrix).getSolver() 212 .solve(b2.matrix); 213 return CommonsMathDenseDoubleMatrix2DFactory.INSTANCE 214 .dense(ret); 215 } else { 216 RealMatrix ret = new QRDecompositionImpl(matrix).getSolver() 217 .solve(b2.matrix); 218 return CommonsMathDenseDoubleMatrix2DFactory.INSTANCE 219 .dense(ret); 220 } 221 } else { 222 return super.solve(b); 223 } 224 } 225 226 public Matrix solveSPD(Matrix b) { 227 try { 228 if (b instanceof AbstractCommonsMathDenseDoubleMatrix2D) { 229 AbstractCommonsMathDenseDoubleMatrix2D b2 = (AbstractCommonsMathDenseDoubleMatrix2D) b; 230 RealMatrix ret = new CholeskyDecompositionImpl(matrix) 231 .getSolver().solve(b2.matrix); 232 return CommonsMathDenseDoubleMatrix2DFactory.INSTANCE 233 .dense(ret); 234 } else { 235 return super.solve(b); 236 } 237 } catch (Exception e) { 238 throw new MatrixException(e); 239 } 240 } 241 }