Merge pull request #189 from mipt-npm/commandertvis/cm-decompositions

Support of new matrix decomposition features for CMMatrix

kmath-commons/src/main/kotlin/kscience/kmath/commons/linear/CMMatrix.kt

@@ -1,11 +1,9 @@
 package kscience.kmath.commons.linear
 
-import kscience.kmath.linear.DiagonalFeature
+import kscience.kmath.linear.*
-import kscience.kmath.linear.MatrixContext
-import kscience.kmath.linear.Point
-import kscience.kmath.linear.origin
 import kscience.kmath.misc.UnstableKMathAPI
 import kscience.kmath.structures.Matrix
+import kscience.kmath.structures.RealBuffer
 import org.apache.commons.math3.linear.*
 import kotlin.reflect.KClass
 import kotlin.reflect.cast
@@ -17,8 +15,40 @@ public inline class CMMatrix(public val origin: RealMatrix) : Matrix<Double> {
     @UnstableKMathAPI
     override fun <T : Any> getFeature(type: KClass<T>): T? = when (type) {
         DiagonalFeature::class -> if (origin is DiagonalMatrix) DiagonalFeature else null
+
+        DeterminantFeature::class, LupDecompositionFeature::class -> object :
+            DeterminantFeature<Double>,
+            LupDecompositionFeature<Double> {
+            private val lup by lazy { LUDecomposition(origin) }
+            override val determinant: Double by lazy { lup.determinant }
+            override val l: Matrix<Double> by lazy { CMMatrix(lup.l) + LFeature }
+            override val u: Matrix<Double> by lazy { CMMatrix(lup.u) + UFeature }
+            override val p: Matrix<Double> by lazy { CMMatrix(lup.p) }
+        }
+
+        CholeskyDecompositionFeature::class -> object : CholeskyDecompositionFeature<Double> {
+            override val l: Matrix<Double> by lazy {
+                val cholesky = CholeskyDecomposition(origin)
+                CMMatrix(cholesky.l) + LFeature
+            }
+        }
+
+        QRDecompositionFeature::class -> object : QRDecompositionFeature<Double> {
+            private val qr by lazy { QRDecomposition(origin) }
+            override val q: Matrix<Double> by lazy { CMMatrix(qr.q) + OrthogonalFeature }
+            override val r: Matrix<Double> by lazy { CMMatrix(qr.r) + UFeature }
+        }
+
+        SingularValueDecompositionFeature::class -> object : SingularValueDecompositionFeature<Double> {
+            private val sv by lazy { SingularValueDecomposition(origin) }
+            override val u: Matrix<Double> by lazy { CMMatrix(sv.u) }
+            override val s: Matrix<Double> by lazy { CMMatrix(sv.s) }
+            override val v: Matrix<Double> by lazy { CMMatrix(sv.v) }
+            override val singularValues: Point<Double> by lazy { RealBuffer(sv.singularValues) }
+        }
+
         else -> null
-    }?.let { type.cast(it) }
+    }?.let(type::cast)
 
     public override operator fun get(i: Int, j: Int): Double = origin.getEntry(i, j)
 }

kmath-core/src/commonMain/kotlin/kscience/kmath/linear/MatrixFeatures.kt

@@ -39,6 +39,8 @@ public interface InverseMatrixFeature<T : Any> : MatrixFeature {
 
 /**
  * Matrices with this feature can compute their determinant.
+ *
+ * @param T the type of matrices' items.
  */
 public interface DeterminantFeature<T : Any> : MatrixFeature {
     /**

kmath-ejml/src/main/kotlin/kscience/kmath/ejml/EjmlMatrix.kt

@@ -16,21 +16,20 @@ import kotlin.reflect.cast
  * @property origin the underlying [SimpleMatrix].
  * @author Iaroslav Postovalov
  */
-public class EjmlMatrix(
+public class EjmlMatrix(public val origin: SimpleMatrix) : Matrix<Double> {
-    public val origin: SimpleMatrix,
-) : Matrix<Double> {
     public override val rowNum: Int get() = origin.numRows()
-
     public override val colNum: Int get() = origin.numCols()
 
     @UnstableKMathAPI
-    override fun <T : Any> getFeature(type: KClass<T>): T? = when (type) {
+    public override fun <T : Any> getFeature(type: KClass<T>): T? = when (type) {
         InverseMatrixFeature::class -> object : InverseMatrixFeature<Double> {
             override val inverse: Matrix<Double> by lazy { EjmlMatrix(origin.invert()) }
         }
+
         DeterminantFeature::class -> object : DeterminantFeature<Double> {
             override val determinant: Double by lazy(origin::determinant)
         }
+
         SingularValueDecompositionFeature::class -> object : SingularValueDecompositionFeature<Double> {
             private val svd by lazy {
                 DecompositionFactory_DDRM.svd(origin.numRows(), origin.numCols(), true, true, false)
@@ -42,14 +41,19 @@ public class EjmlMatrix(
             override val v: Matrix<Double> by lazy { EjmlMatrix(SimpleMatrix(svd.getV(null, false))) }
             override val singularValues: Point<Double> by lazy { RealBuffer(svd.singularValues) }
         }
+
         QRDecompositionFeature::class -> object : QRDecompositionFeature<Double> {
             private val qr by lazy {
                 DecompositionFactory_DDRM.qr().apply { decompose(origin.ddrm.copy()) }
             }
 
-            override val q: Matrix<Double> by lazy { EjmlMatrix(SimpleMatrix(qr.getQ(null, false))) }
+            override val q: Matrix<Double> by lazy {
-            override val r: Matrix<Double> by lazy { EjmlMatrix(SimpleMatrix(qr.getR(null, false))) }
+                EjmlMatrix(SimpleMatrix(qr.getQ(null, false))) + OrthogonalFeature
             }
+
+            override val r: Matrix<Double> by lazy { EjmlMatrix(SimpleMatrix(qr.getR(null, false))) + UFeature }
+        }
+
         CholeskyDecompositionFeature::class -> object : CholeskyDecompositionFeature<Double> {
             override val l: Matrix<Double> by lazy {
                 val cholesky =
@@ -58,6 +62,7 @@ public class EjmlMatrix(
                 EjmlMatrix(SimpleMatrix(cholesky.getT(null))) + LFeature
             }
         }
+
         LupDecompositionFeature::class -> object : LupDecompositionFeature<Double> {
             private val lup by lazy {
                 DecompositionFactory_DDRM.lu(origin.numRows(), origin.numCols()).apply { decompose(origin.ddrm.copy()) }
@@ -73,8 +78,9 @@ public class EjmlMatrix(
 
             override val p: Matrix<Double> by lazy { EjmlMatrix(SimpleMatrix(lup.getRowPivot(null))) }
         }
+
         else -> null
-    }?.let { type.cast(it) }
+    }?.let(type::cast)
 
     public override operator fun get(i: Int, j: Int): Double = origin[i, j]