com.opengamma.strata.math.impl.matrix.CommonsMatrixAlgebra Maven / Gradle / Ivy

Go to download
Show more of this group Show more artifacts with this name
Show all versions of strata-math Show documentation
Mathematic support for Strata
There is a newer version: 2.12.46
/*
 * Copyright (C) 2009 - present by OpenGamma Inc. and the OpenGamma group of companies
 *
 * Please see distribution for license.
 */
package com.opengamma.strata.math.impl.matrix;

import org.apache.commons.math3.linear.LUDecomposition;
import org.apache.commons.math3.linear.RealMatrix;
import org.apache.commons.math3.linear.RealVector;
import org.apache.commons.math3.linear.SingularValueDecomposition;

import com.opengamma.strata.collect.ArgChecker;
import com.opengamma.strata.collect.array.DoubleArray;
import com.opengamma.strata.collect.array.DoubleMatrix;
import com.opengamma.strata.collect.array.Matrix;
import com.opengamma.strata.math.impl.util.CommonsMathWrapper;

/**
 * Provides matrix algebra by using the Commons library. 
 */
public class CommonsMatrixAlgebra extends MatrixAlgebra {

  @Override
  public double getCondition(Matrix m) {
    ArgChecker.notNull(m, "m");
    if (m instanceof DoubleMatrix) {
      RealMatrix temp = CommonsMathWrapper.wrap((DoubleMatrix) m);
      SingularValueDecomposition svd = new SingularValueDecomposition(temp);
      return svd.getConditionNumber();
    }
    throw new IllegalArgumentException("Can only find condition number of DoubleMatrix; have " + m.getClass());
  }

  @Override
  public double getDeterminant(Matrix m) {
    ArgChecker.notNull(m, "m");
    if (m instanceof DoubleMatrix) {
      RealMatrix temp = CommonsMathWrapper.wrap((DoubleMatrix) m);
      LUDecomposition lud = new LUDecomposition(temp);
      return lud.getDeterminant();
    }
    throw new IllegalArgumentException("Can only find determinant of DoubleMatrix; have " + m.getClass());
  }

  @Override
  public double getInnerProduct(Matrix m1, Matrix m2) {
    ArgChecker.notNull(m1, "m1");
    ArgChecker.notNull(m2, "m2");
    if (m1 instanceof DoubleArray && m2 instanceof DoubleArray) {
      RealVector t1 = CommonsMathWrapper.wrap((DoubleArray) m1);
      RealVector t2 = CommonsMathWrapper.wrap((DoubleArray) m2);
      return t1.dotProduct(t2);
    }
    throw new IllegalArgumentException("Can only find inner product of DoubleArray; have " +
        m1.getClass() + " and " + m2.getClass());
  }

  @Override
  public DoubleMatrix getInverse(Matrix m) {
    ArgChecker.notNull(m, "matrix was null");
    if (m instanceof DoubleMatrix) {
      RealMatrix temp = CommonsMathWrapper.wrap((DoubleMatrix) m);
      SingularValueDecomposition sv = new SingularValueDecomposition(temp);
      RealMatrix inv = sv.getSolver().getInverse();
      return CommonsMathWrapper.unwrap(inv);
    }
    throw new IllegalArgumentException("Can only find inverse of DoubleMatrix; have " + m.getClass());
  }

  @Override
  public double getNorm1(Matrix m) {
    ArgChecker.notNull(m, "m");
    if (m instanceof DoubleArray) {
      RealVector temp = CommonsMathWrapper.wrap((DoubleArray) m);
      return temp.getL1Norm();
    } else if (m instanceof DoubleMatrix) {
      RealMatrix temp = CommonsMathWrapper.wrap((DoubleMatrix) m);
      // TODO find if commons implements this anywhere, so we are not doing it
      // by hand
      double max = 0.0;
      for (int col = temp.getColumnDimension() - 1; col >= 0; col--) {
        max = Math.max(max, temp.getColumnVector(col).getL1Norm());
      }
      return max;

    }
    throw new IllegalArgumentException("Can only find norm1 of DoubleMatrix; have " + m.getClass());
  }

  @Override
  public double getNorm2(Matrix m) {
    ArgChecker.notNull(m, "m");
    if (m instanceof DoubleArray) {
      RealVector temp = CommonsMathWrapper.wrap((DoubleArray) m);
      return temp.getNorm();
    } else if (m instanceof DoubleMatrix) {
      RealMatrix temp = CommonsMathWrapper.wrap((DoubleMatrix) m);
      SingularValueDecomposition svd = new SingularValueDecomposition(temp);
      return svd.getNorm();
    }
    throw new IllegalArgumentException("Can only find norm2 of DoubleMatrix; have " + m.getClass());
  }

  @Override
  public double getNormInfinity(Matrix m) {
    ArgChecker.notNull(m, "m");
    if (m instanceof DoubleArray) {
      RealVector temp = CommonsMathWrapper.wrap((DoubleArray) m);
      return temp.getLInfNorm();
    } else if (m instanceof DoubleMatrix) {
      RealMatrix temp = CommonsMathWrapper.wrap((DoubleMatrix) m);
      //REVIEW Commons getNorm() is wrong - it returns the column norm
      // TODO find if commons implements this anywhere, so we are not doing it
      // by hand
      double max = 0.0;
      for (int row = temp.getRowDimension() - 1; row >= 0; row--) {
        max = Math.max(max, temp.getRowVector(row).getL1Norm());
      }
      return max;
    }
    throw new IllegalArgumentException("Can only find normInfinity of DoubleMatrix; have " + m.getClass());
  }

  @Override
  public DoubleMatrix getOuterProduct(Matrix m1, Matrix m2) {
    ArgChecker.notNull(m1, "m1");
    ArgChecker.notNull(m2, "m2");
    if (m1 instanceof DoubleArray && m2 instanceof DoubleArray) {
      RealVector t1 = CommonsMathWrapper.wrap((DoubleArray) m1);
      RealVector t2 = CommonsMathWrapper.wrap((DoubleArray) m2);
      return CommonsMathWrapper.unwrap(t1.outerProduct(t2));
    }
    throw new IllegalArgumentException("Can only find outer product of DoubleArray; have " +
        m1.getClass() + " and " + m2.getClass());
  }

  @Override
  public DoubleMatrix getPower(Matrix m, int p) {
    ArgChecker.notNull(m, "m");
    RealMatrix temp;
    if (m instanceof DoubleMatrix) {
      temp = CommonsMathWrapper.wrap((DoubleMatrix) m);
    } else {
      throw new IllegalArgumentException("Can only find powers of DoubleMatrix; have " + m.getClass());
    }
    return CommonsMathWrapper.unwrap(temp.power(p));
  }

  /**
   * Returns a real matrix raised to some real power 
   * Currently this method is limited to symmetric matrices only as Commons Math does not
   * support the diagonalization of asymmetric matrices.
   * 
   * @param m The symmetric matrix to take the power of. 
   * @param p The power to raise to matrix to
   * @return The result
   */
  @Override
  public DoubleMatrix getPower(Matrix m, double p) {
    throw new UnsupportedOperationException();
  }

  @Override
  public double getTrace(Matrix m) {
    ArgChecker.notNull(m, "m");
    if (m instanceof DoubleMatrix) {
      RealMatrix temp = CommonsMathWrapper.wrap((DoubleMatrix) m);
      return temp.getTrace();
    }
    throw new IllegalArgumentException("Can only find trace of DoubleMatrix; have " + m.getClass());
  }

  @Override
  public DoubleMatrix getTranspose(Matrix m) {
    ArgChecker.notNull(m, "m");
    if (m instanceof DoubleMatrix) {
      RealMatrix temp = CommonsMathWrapper.wrap((DoubleMatrix) m);
      return CommonsMathWrapper.unwrap(temp.transpose());
    }
    throw new IllegalArgumentException("Can only find transpose of DoubleMatrix; have " + m.getClass());
  }

  /**
   * {@inheritDoc}
   * The following combinations of input matrices m1 and m2 are allowed:
   * 
   *  m1 = 2-D matrix, m2 = 2-D matrix, returns $\mathbf{C} = \mathbf{AB}$
   * 
 m1 = 2-D matrix, m2 = 1-D matrix, returns $\mathbf{C} = \mathbf{A}b$
   * 
   */
  @Override
  public Matrix multiply(Matrix m1, Matrix m2) {
    ArgChecker.notNull(m1, "m1");
    ArgChecker.notNull(m2, "m2");
    ArgChecker.isTrue(!(m1 instanceof DoubleArray), "Cannot have 1D matrix as first argument");
    if (m1 instanceof DoubleMatrix) {
      RealMatrix t1 = CommonsMathWrapper.wrap((DoubleMatrix) m1);
      RealMatrix t2;
      if (m2 instanceof DoubleArray) {
        t2 = CommonsMathWrapper.wrapAsMatrix((DoubleArray) m2);
      } else if (m2 instanceof DoubleMatrix) {
        t2 = CommonsMathWrapper.wrap((DoubleMatrix) m2);
      } else {
        throw new IllegalArgumentException("Can only have 1D or 2D matrix as second argument");
      }
      return CommonsMathWrapper.unwrap(t1.multiply(t2));
    }
    throw new IllegalArgumentException("Can only multiply 2D and 1D matrices");
  }

}