forked from kscience/kmath
Merge pull request #339 from mipt-npm/commandertvis/ejml
Code generation of EJML linear spaces
This commit is contained in:
Iaroslav Postovalov
GitHub
commit
47526d4dd5
5
buildSrc/build.gradle.kts
Normal file
5
buildSrc/build.gradle.kts
Normal file
@ -0,0 +1,5 @@
|
||||
plugins {
|
||||
`kotlin-dsl`
|
||||
}
|
||||
|
||||
repositories.mavenCentral()
|
||||
@ -0,0 +1,415 @@
|
||||
/*
|
||||
* Copyright 2018-2021 KMath contributors.
|
||||
* Use of this source code is governed by the Apache 2.0 license that can be found in the license/LICENSE.txt file.
|
||||
*/
|
||||
|
||||
@file:Suppress("KDocUnresolvedReference")
|
||||
|
||||
package space.kscience.kmath.ejml.codegen
|
||||
|
||||
import org.intellij.lang.annotations.Language
|
||||
import java.io.File
|
||||
|
||||
private fun Appendable.appendEjmlVector(type: String, ejmlMatrixType: String) {
|
||||
@Language("kotlin") val text = """/**
|
||||
* [EjmlVector] specialization for [$type].
|
||||
*/
|
||||
public class Ejml${type}Vector<out M : $ejmlMatrixType>(public override val origin: M) : EjmlVector<$type, M>(origin) {
|
||||
init {
|
||||
require(origin.numRows == 1) { "The origin matrix must have only one row to form a vector" }
|
||||
}
|
||||
|
||||
public override operator fun get(index: Int): $type = origin[0, index]
|
||||
}"""
|
||||
appendLine(text)
|
||||
appendLine()
|
||||
}
|
||||
|
||||
private fun Appendable.appendEjmlMatrix(type: String, ejmlMatrixType: String) {
|
||||
val text = """/**
|
||||
* [EjmlMatrix] specialization for [$type].
|
||||
*/
|
||||
public class Ejml${type}Matrix<out M : $ejmlMatrixType>(public override val origin: M) : EjmlMatrix<$type, M>(origin) {
|
||||
public override operator fun get(i: Int, j: Int): $type = origin[i, j]
|
||||
}"""
|
||||
appendLine(text)
|
||||
appendLine()
|
||||
}
|
||||
|
||||
private fun Appendable.appendEjmlLinearSpace(
|
||||
type: String,
|
||||
kmathAlgebra: String,
|
||||
ejmlMatrixParentTypeMatrix: String,
|
||||
ejmlMatrixType: String,
|
||||
ejmlMatrixDenseType: String,
|
||||
ops: String,
|
||||
denseOps: String,
|
||||
isDense: Boolean,
|
||||
) {
|
||||
@Language("kotlin") val text = """/**
|
||||
* [EjmlLinearSpace] implementation based on [CommonOps_$ops], [DecompositionFactory_${ops}] operations and
|
||||
* [${ejmlMatrixType}] matrices.
|
||||
*/
|
||||
public object EjmlLinearSpace${ops} : EjmlLinearSpace<${type}, ${kmathAlgebra}, $ejmlMatrixType>() {
|
||||
/**
|
||||
* The [${kmathAlgebra}] reference.
|
||||
*/
|
||||
public override val elementAlgebra: $kmathAlgebra get() = $kmathAlgebra
|
||||
|
||||
@Suppress("UNCHECKED_CAST")
|
||||
public override fun Matrix<${type}>.toEjml(): Ejml${type}Matrix<${ejmlMatrixType}> = when {
|
||||
this is Ejml${type}Matrix<*> && origin is $ejmlMatrixType -> this as Ejml${type}Matrix<${ejmlMatrixType}>
|
||||
else -> buildMatrix(rowNum, colNum) { i, j -> get(i, j) }
|
||||
}
|
||||
|
||||
@Suppress("UNCHECKED_CAST")
|
||||
public override fun Point<${type}>.toEjml(): Ejml${type}Vector<${ejmlMatrixType}> = when {
|
||||
this is Ejml${type}Vector<*> && origin is $ejmlMatrixType -> this as Ejml${type}Vector<${ejmlMatrixType}>
|
||||
else -> Ejml${type}Vector(${ejmlMatrixType}(size, 1).also {
|
||||
(0 until it.numRows).forEach { row -> it[row, 0] = get(row) }
|
||||
})
|
||||
}
|
||||
|
||||
public override fun buildMatrix(
|
||||
rows: Int,
|
||||
columns: Int,
|
||||
initializer: ${kmathAlgebra}.(i: Int, j: Int) -> ${type},
|
||||
): Ejml${type}Matrix<${ejmlMatrixType}> = ${ejmlMatrixType}(rows, columns).also {
|
||||
(0 until rows).forEach { row ->
|
||||
(0 until columns).forEach { col -> it[row, col] = elementAlgebra.initializer(row, col) }
|
||||
}
|
||||
}.wrapMatrix()
|
||||
|
||||
public override fun buildVector(
|
||||
size: Int,
|
||||
initializer: ${kmathAlgebra}.(Int) -> ${type},
|
||||
): Ejml${type}Vector<${ejmlMatrixType}> = Ejml${type}Vector(${ejmlMatrixType}(size, 1).also {
|
||||
(0 until it.numRows).forEach { row -> it[row, 0] = elementAlgebra.initializer(row) }
|
||||
})
|
||||
|
||||
private fun <T : ${ejmlMatrixParentTypeMatrix}> T.wrapMatrix() = Ejml${type}Matrix(this)
|
||||
private fun <T : ${ejmlMatrixParentTypeMatrix}> T.wrapVector() = Ejml${type}Vector(this)
|
||||
|
||||
public override fun Matrix<${type}>.unaryMinus(): Matrix<${type}> = this * elementAlgebra { -one }
|
||||
|
||||
public override fun Matrix<${type}>.dot(other: Matrix<${type}>): Ejml${type}Matrix<${ejmlMatrixType}> {
|
||||
val out = ${ejmlMatrixType}(1, 1)
|
||||
CommonOps_${ops}.mult(toEjml().origin, other.toEjml().origin, out)
|
||||
return out.wrapMatrix()
|
||||
}
|
||||
|
||||
public override fun Matrix<${type}>.dot(vector: Point<${type}>): Ejml${type}Vector<${ejmlMatrixType}> {
|
||||
val out = ${ejmlMatrixType}(1, 1)
|
||||
CommonOps_${ops}.mult(toEjml().origin, vector.toEjml().origin, out)
|
||||
return out.wrapVector()
|
||||
}
|
||||
|
||||
public override operator fun Matrix<${type}>.minus(other: Matrix<${type}>): Ejml${type}Matrix<${ejmlMatrixType}> {
|
||||
val out = ${ejmlMatrixType}(1, 1)
|
||||
|
||||
CommonOps_${ops}.add(
|
||||
elementAlgebra.one,
|
||||
toEjml().origin,
|
||||
elementAlgebra { -one },
|
||||
other.toEjml().origin,
|
||||
out,${
|
||||
if (isDense) "" else
|
||||
"""
|
||||
null,
|
||||
null,"""
|
||||
}
|
||||
)
|
||||
|
||||
return out.wrapMatrix()
|
||||
}
|
||||
|
||||
public override operator fun Matrix<${type}>.times(value: ${type}): Ejml${type}Matrix<${ejmlMatrixType}> {
|
||||
val res = ${ejmlMatrixType}(1, 1)
|
||||
CommonOps_${ops}.scale(value, toEjml().origin, res)
|
||||
return res.wrapMatrix()
|
||||
}
|
||||
|
||||
public override fun Point<${type}>.unaryMinus(): Ejml${type}Vector<${ejmlMatrixType}> {
|
||||
val res = ${ejmlMatrixType}(1, 1)
|
||||
CommonOps_${ops}.changeSign(toEjml().origin, res)
|
||||
return res.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Matrix<${type}>.plus(other: Matrix<${type}>): Ejml${type}Matrix<${ejmlMatrixType}> {
|
||||
val out = ${ejmlMatrixType}(1, 1)
|
||||
|
||||
CommonOps_${ops}.add(
|
||||
elementAlgebra.one,
|
||||
toEjml().origin,
|
||||
elementAlgebra.one,
|
||||
other.toEjml().origin,
|
||||
out,${
|
||||
if (isDense) "" else
|
||||
"""
|
||||
null,
|
||||
null,"""
|
||||
}
|
||||
)
|
||||
|
||||
return out.wrapMatrix()
|
||||
}
|
||||
|
||||
public override fun Point<${type}>.plus(other: Point<${type}>): Ejml${type}Vector<${ejmlMatrixType}> {
|
||||
val out = ${ejmlMatrixType}(1, 1)
|
||||
|
||||
CommonOps_${ops}.add(
|
||||
elementAlgebra.one,
|
||||
toEjml().origin,
|
||||
elementAlgebra.one,
|
||||
other.toEjml().origin,
|
||||
out,${
|
||||
if (isDense) "" else
|
||||
"""
|
||||
null,
|
||||
null,"""
|
||||
}
|
||||
)
|
||||
|
||||
return out.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Point<${type}>.minus(other: Point<${type}>): Ejml${type}Vector<${ejmlMatrixType}> {
|
||||
val out = ${ejmlMatrixType}(1, 1)
|
||||
|
||||
CommonOps_${ops}.add(
|
||||
elementAlgebra.one,
|
||||
toEjml().origin,
|
||||
elementAlgebra { -one },
|
||||
other.toEjml().origin,
|
||||
out,${
|
||||
if (isDense) "" else
|
||||
"""
|
||||
null,
|
||||
null,"""
|
||||
}
|
||||
)
|
||||
|
||||
return out.wrapVector()
|
||||
}
|
||||
|
||||
public override fun ${type}.times(m: Matrix<${type}>): Ejml${type}Matrix<${ejmlMatrixType}> = m * this
|
||||
|
||||
public override fun Point<${type}>.times(value: ${type}): Ejml${type}Vector<${ejmlMatrixType}> {
|
||||
val res = ${ejmlMatrixType}(1, 1)
|
||||
CommonOps_${ops}.scale(value, toEjml().origin, res)
|
||||
return res.wrapVector()
|
||||
}
|
||||
|
||||
public override fun ${type}.times(v: Point<${type}>): Ejml${type}Vector<${ejmlMatrixType}> = v * this
|
||||
|
||||
@UnstableKMathAPI
|
||||
public override fun <F : StructureFeature> getFeature(structure: Matrix<${type}>, type: KClass<out F>): F? {
|
||||
structure.getFeature(type)?.let { return it }
|
||||
val origin = structure.toEjml().origin
|
||||
|
||||
return when (type) {
|
||||
${
|
||||
if (isDense)
|
||||
""" InverseMatrixFeature::class -> object : InverseMatrixFeature<${type}> {
|
||||
override val inverse: Matrix<${type}> by lazy {
|
||||
val res = origin.copy()
|
||||
CommonOps_${ops}.invert(res)
|
||||
res.wrapMatrix()
|
||||
}
|
||||
}
|
||||
|
||||
DeterminantFeature::class -> object : DeterminantFeature<${type}> {
|
||||
override val determinant: $type by lazy { CommonOps_${ops}.det(origin) }
|
||||
}
|
||||
|
||||
SingularValueDecompositionFeature::class -> object : SingularValueDecompositionFeature<${type}> {
|
||||
private val svd by lazy {
|
||||
DecompositionFactory_${ops}.svd(origin.numRows, origin.numCols, true, true, false)
|
||||
.apply { decompose(origin.copy()) }
|
||||
}
|
||||
|
||||
override val u: Matrix<${type}> by lazy { svd.getU(null, false).wrapMatrix() }
|
||||
override val s: Matrix<${type}> by lazy { svd.getW(null).wrapMatrix() }
|
||||
override val v: Matrix<${type}> by lazy { svd.getV(null, false).wrapMatrix() }
|
||||
override val singularValues: Point<${type}> by lazy { ${type}Buffer(svd.singularValues) }
|
||||
}
|
||||
|
||||
QRDecompositionFeature::class -> object : QRDecompositionFeature<${type}> {
|
||||
private val qr by lazy {
|
||||
DecompositionFactory_${ops}.qr().apply { decompose(origin.copy()) }
|
||||
}
|
||||
|
||||
override val q: Matrix<${type}> by lazy {
|
||||
qr.getQ(null, false).wrapMatrix() + OrthogonalFeature
|
||||
}
|
||||
|
||||
override val r: Matrix<${type}> by lazy { qr.getR(null, false).wrapMatrix() + UFeature }
|
||||
}
|
||||
|
||||
CholeskyDecompositionFeature::class -> object : CholeskyDecompositionFeature<${type}> {
|
||||
override val l: Matrix<${type}> by lazy {
|
||||
val cholesky =
|
||||
DecompositionFactory_${ops}.chol(structure.rowNum, true).apply { decompose(origin.copy()) }
|
||||
|
||||
cholesky.getT(null).wrapMatrix() + LFeature
|
||||
}
|
||||
}
|
||||
|
||||
LupDecompositionFeature::class -> object : LupDecompositionFeature<${type}> {
|
||||
private val lup by lazy {
|
||||
DecompositionFactory_${ops}.lu(origin.numRows, origin.numCols).apply { decompose(origin.copy()) }
|
||||
}
|
||||
|
||||
override val l: Matrix<${type}> by lazy {
|
||||
lup.getLower(null).wrapMatrix() + LFeature
|
||||
}
|
||||
|
||||
override val u: Matrix<${type}> by lazy {
|
||||
lup.getUpper(null).wrapMatrix() + UFeature
|
||||
}
|
||||
|
||||
override val p: Matrix<${type}> by lazy { lup.getRowPivot(null).wrapMatrix() }
|
||||
}""" else """ QRDecompositionFeature::class -> object : QRDecompositionFeature<$type> {
|
||||
private val qr by lazy {
|
||||
DecompositionFactory_${ops}.qr(FillReducing.NONE).apply { decompose(origin.copy()) }
|
||||
}
|
||||
|
||||
override val q: Matrix<${type}> by lazy {
|
||||
qr.getQ(null, false).wrapMatrix() + OrthogonalFeature
|
||||
}
|
||||
|
||||
override val r: Matrix<${type}> by lazy { qr.getR(null, false).wrapMatrix() + UFeature }
|
||||
}
|
||||
|
||||
CholeskyDecompositionFeature::class -> object : CholeskyDecompositionFeature<${type}> {
|
||||
override val l: Matrix<${type}> by lazy {
|
||||
val cholesky =
|
||||
DecompositionFactory_${ops}.cholesky().apply { decompose(origin.copy()) }
|
||||
|
||||
(cholesky.getT(null) as ${ejmlMatrixParentTypeMatrix}).wrapMatrix() + LFeature
|
||||
}
|
||||
}
|
||||
|
||||
LUDecompositionFeature::class, DeterminantFeature::class, InverseMatrixFeature::class -> object :
|
||||
LUDecompositionFeature<${type}>, DeterminantFeature<${type}>, InverseMatrixFeature<${type}> {
|
||||
private val lu by lazy {
|
||||
DecompositionFactory_${ops}.lu(FillReducing.NONE).apply { decompose(origin.copy()) }
|
||||
}
|
||||
|
||||
override val l: Matrix<${type}> by lazy {
|
||||
lu.getLower(null).wrapMatrix() + LFeature
|
||||
}
|
||||
|
||||
override val u: Matrix<${type}> by lazy {
|
||||
lu.getUpper(null).wrapMatrix() + UFeature
|
||||
}
|
||||
|
||||
override val inverse: Matrix<${type}> by lazy {
|
||||
var a = origin
|
||||
val inverse = ${ejmlMatrixDenseType}(1, 1)
|
||||
val solver = LinearSolverFactory_${ops}.lu(FillReducing.NONE)
|
||||
if (solver.modifiesA()) a = a.copy()
|
||||
val i = CommonOps_${denseOps}.identity(a.numRows)
|
||||
solver.solve(i, inverse)
|
||||
inverse.wrapMatrix()
|
||||
}
|
||||
|
||||
override val determinant: $type by lazy { elementAlgebra.number(lu.computeDeterminant().real) }
|
||||
}"""
|
||||
}
|
||||
|
||||
else -> null
|
||||
}?.let(type::cast)
|
||||
}
|
||||
|
||||
/**
|
||||
* Solves for *x* in the following equation: *x = [a] <sup>-1</sup> · [b]*.
|
||||
*
|
||||
* @param a the base matrix.
|
||||
* @param b n by p matrix.
|
||||
* @return the solution for *x* that is n by p.
|
||||
*/
|
||||
public fun solve(a: Matrix<${type}>, b: Matrix<${type}>): Ejml${type}Matrix<${ejmlMatrixType}> {
|
||||
val res = ${ejmlMatrixType}(1, 1)
|
||||
CommonOps_${ops}.solve(${ejmlMatrixType}(a.toEjml().origin), ${ejmlMatrixType}(b.toEjml().origin), res)
|
||||
return res.wrapMatrix()
|
||||
}
|
||||
|
||||
/**
|
||||
* Solves for *x* in the following equation: *x = [a] <sup>-1</sup> · [b]*.
|
||||
*
|
||||
* @param a the base matrix.
|
||||
* @param b n by p vector.
|
||||
* @return the solution for *x* that is n by p.
|
||||
*/
|
||||
public fun solve(a: Matrix<${type}>, b: Point<${type}>): Ejml${type}Vector<${ejmlMatrixType}> {
|
||||
val res = ${ejmlMatrixType}(1, 1)
|
||||
CommonOps_${ops}.solve(${ejmlMatrixType}(a.toEjml().origin), ${ejmlMatrixType}(b.toEjml().origin), res)
|
||||
return Ejml${type}Vector(res)
|
||||
}
|
||||
}"""
|
||||
appendLine(text)
|
||||
appendLine()
|
||||
}
|
||||
|
||||
|
||||
/**
|
||||
* Generates routine EJML classes.
|
||||
*/
|
||||
fun ejmlCodegen(outputFile: String): Unit = File(outputFile).run {
|
||||
parentFile.mkdirs()
|
||||
|
||||
writer().use {
|
||||
it.appendLine("/*")
|
||||
it.appendLine(" * Copyright 2018-2021 KMath contributors.")
|
||||
it.appendLine(" * Use of this source code is governed by the Apache 2.0 license that can be found in the license/LICENSE.txt file.")
|
||||
it.appendLine(" */")
|
||||
it.appendLine()
|
||||
it.appendLine("/* This file is generated with buildSrc/src/main/kotlin/space/kscience/kmath/ejml/codegen/ejmlCodegen.kt */")
|
||||
it.appendLine()
|
||||
it.appendLine("package space.kscience.kmath.ejml")
|
||||
it.appendLine()
|
||||
it.appendLine("import org.ejml.data.*")
|
||||
it.appendLine("import space.kscience.kmath.linear.*")
|
||||
it.appendLine("import space.kscience.kmath.operations.*")
|
||||
it.appendLine("import space.kscience.kmath.structures.*")
|
||||
it.appendLine("import space.kscience.kmath.misc.*")
|
||||
it.appendLine("import kotlin.reflect.*")
|
||||
it.appendLine("import org.ejml.dense.row.*")
|
||||
it.appendLine("import org.ejml.dense.row.factory.*")
|
||||
it.appendLine("import org.ejml.sparse.*")
|
||||
it.appendLine("import org.ejml.sparse.csc.*")
|
||||
it.appendLine("import org.ejml.sparse.csc.factory.*")
|
||||
it.appendLine("import space.kscience.kmath.nd.*")
|
||||
it.appendLine("import space.kscience.kmath.linear.Matrix")
|
||||
it.appendLine()
|
||||
it.appendEjmlVector("Double", "DMatrix")
|
||||
it.appendEjmlVector("Float", "FMatrix")
|
||||
it.appendEjmlMatrix("Double", "DMatrix")
|
||||
it.appendEjmlMatrix("Float", "FMatrix")
|
||||
it.appendEjmlLinearSpace("Double", "DoubleField", "DMatrix", "DMatrixRMaj", "DMatrixRMaj", "DDRM", "DDRM", true)
|
||||
it.appendEjmlLinearSpace("Float", "FloatField", "FMatrix", "FMatrixRMaj", "FMatrixRMaj", "FDRM", "FDRM", true)
|
||||
|
||||
it.appendEjmlLinearSpace(
|
||||
type = "Double",
|
||||
kmathAlgebra = "DoubleField",
|
||||
ejmlMatrixParentTypeMatrix = "DMatrix",
|
||||
ejmlMatrixType = "DMatrixSparseCSC",
|
||||
ejmlMatrixDenseType = "DMatrixRMaj",
|
||||
ops = "DSCC",
|
||||
denseOps = "DDRM",
|
||||
isDense = false,
|
||||
)
|
||||
|
||||
it.appendEjmlLinearSpace(
|
||||
type = "Float",
|
||||
kmathAlgebra = "FloatField",
|
||||
ejmlMatrixParentTypeMatrix = "FMatrix",
|
||||
ejmlMatrixType = "FMatrixSparseCSC",
|
||||
ejmlMatrixDenseType = "FMatrixRMaj",
|
||||
ops = "FSCC",
|
||||
denseOps = "FDRM",
|
||||
isDense = false,
|
||||
)
|
||||
}
|
||||
}
|
||||
@ -75,6 +75,23 @@ public object LFeature : MatrixFeature
|
||||
*/
|
||||
public object UFeature : MatrixFeature
|
||||
|
||||
/**
|
||||
* Matrices with this feature support LU factorization: *a = [l] · [u]* where *a* is the owning matrix.
|
||||
*
|
||||
* @param T the type of matrices' items.
|
||||
*/
|
||||
public interface LUDecompositionFeature<T : Any> : MatrixFeature {
|
||||
/**
|
||||
* The lower triangular matrix in this decomposition. It may have [LFeature].
|
||||
*/
|
||||
public val l: Matrix<T>
|
||||
|
||||
/**
|
||||
* The upper triangular matrix in this decomposition. It may have [UFeature].
|
||||
*/
|
||||
public val u: Matrix<T>
|
||||
}
|
||||
|
||||
/**
|
||||
* Matrices with this feature support LU factorization with partial pivoting: *[p] · a = [l] · [u]* where
|
||||
* *a* is the owning matrix.
|
||||
|
||||
@ -1,3 +1,5 @@
|
||||
import space.kscience.kmath.ejml.codegen.ejmlCodegen
|
||||
|
||||
plugins {
|
||||
kotlin("jvm")
|
||||
id("ru.mipt.npm.gradle.common")
|
||||
@ -5,6 +7,9 @@ plugins {
|
||||
|
||||
dependencies {
|
||||
api("org.ejml:ejml-ddense:0.40")
|
||||
api("org.ejml:ejml-fdense:0.40")
|
||||
api("org.ejml:ejml-dsparse:0.40")
|
||||
api("org.ejml:ejml-fsparse:0.40")
|
||||
api(project(":kmath-core"))
|
||||
}
|
||||
|
||||
@ -14,19 +19,24 @@ readme {
|
||||
|
||||
feature(
|
||||
id = "ejml-vector",
|
||||
description = "Point implementations.",
|
||||
ref = "src/main/kotlin/space/kscience/kmath/ejml/EjmlVector.kt"
|
||||
)
|
||||
) { "Point implementations." }
|
||||
|
||||
feature(
|
||||
id = "ejml-matrix",
|
||||
description = "Matrix implementation.",
|
||||
ref = "src/main/kotlin/space/kscience/kmath/ejml/EjmlMatrix.kt"
|
||||
)
|
||||
) { "Matrix implementation." }
|
||||
|
||||
feature(
|
||||
id = "ejml-linear-space",
|
||||
description = "LinearSpace implementations.",
|
||||
ref = "src/main/kotlin/space/kscience/kmath/ejml/EjmlLinearSpace.kt"
|
||||
)
|
||||
) { "LinearSpace implementations." }
|
||||
}
|
||||
|
||||
kotlin.sourceSets.main {
|
||||
val codegen by tasks.creating {
|
||||
ejmlCodegen(kotlin.srcDirs.first().absolutePath + "/space/kscience/kmath/ejml/_generated.kt")
|
||||
}
|
||||
|
||||
kotlin.srcDirs(files().builtBy(codegen))
|
||||
}
|
||||
|
||||
@ -5,19 +5,10 @@
|
||||
|
||||
package space.kscience.kmath.ejml
|
||||
|
||||
import org.ejml.data.DMatrix
|
||||
import org.ejml.data.DMatrixD1
|
||||
import org.ejml.data.DMatrixRMaj
|
||||
import org.ejml.dense.row.CommonOps_DDRM
|
||||
import org.ejml.dense.row.factory.DecompositionFactory_DDRM
|
||||
import space.kscience.kmath.linear.*
|
||||
import space.kscience.kmath.misc.UnstableKMathAPI
|
||||
import space.kscience.kmath.nd.StructureFeature
|
||||
import space.kscience.kmath.operations.DoubleField
|
||||
import space.kscience.kmath.linear.LinearSpace
|
||||
import space.kscience.kmath.linear.Matrix
|
||||
import space.kscience.kmath.linear.Point
|
||||
import space.kscience.kmath.operations.Ring
|
||||
import space.kscience.kmath.structures.DoubleBuffer
|
||||
import kotlin.reflect.KClass
|
||||
import kotlin.reflect.cast
|
||||
|
||||
/**
|
||||
* [LinearSpace] implementation specialized for a certain EJML type.
|
||||
@ -27,7 +18,7 @@ import kotlin.reflect.cast
|
||||
* @param M the EJML matrix type.
|
||||
* @author Iaroslav Postovalov
|
||||
*/
|
||||
public abstract class EjmlLinearSpace<T : Any, out A : Ring<T>, M : org.ejml.data.Matrix> : LinearSpace<T, A> {
|
||||
public abstract class EjmlLinearSpace<T : Any, out A : Ring<T>, out M : org.ejml.data.Matrix> : LinearSpace<T, A> {
|
||||
/**
|
||||
* Converts this matrix to EJML one.
|
||||
*/
|
||||
@ -46,209 +37,3 @@ public abstract class EjmlLinearSpace<T : Any, out A : Ring<T>, M : org.ejml.dat
|
||||
|
||||
public abstract override fun buildVector(size: Int, initializer: A.(Int) -> T): EjmlVector<T, M>
|
||||
}
|
||||
|
||||
/**
|
||||
* [EjmlLinearSpace] implementation based on [CommonOps_DDRM], [DecompositionFactory_DDRM] operations and
|
||||
* [DMatrixRMaj] matrices.
|
||||
*
|
||||
* @author Iaroslav Postovalov
|
||||
*/
|
||||
public object EjmlLinearSpaceDDRM : EjmlLinearSpace<Double, DoubleField, DMatrixRMaj>() {
|
||||
/**
|
||||
* The [DoubleField] reference.
|
||||
*/
|
||||
public override val elementAlgebra: DoubleField get() = DoubleField
|
||||
|
||||
@Suppress("UNCHECKED_CAST")
|
||||
public override fun Matrix<Double>.toEjml(): EjmlDoubleMatrix<DMatrixRMaj> = when {
|
||||
this is EjmlDoubleMatrix<*> && origin is DMatrixRMaj -> this as EjmlDoubleMatrix<DMatrixRMaj>
|
||||
else -> buildMatrix(rowNum, colNum) { i, j -> get(i, j) }
|
||||
}
|
||||
|
||||
@Suppress("UNCHECKED_CAST")
|
||||
public override fun Point<Double>.toEjml(): EjmlDoubleVector<DMatrixRMaj> = when {
|
||||
this is EjmlDoubleVector<*> && origin is DMatrixRMaj -> this as EjmlDoubleVector<DMatrixRMaj>
|
||||
else -> EjmlDoubleVector(DMatrixRMaj(size, 1).also {
|
||||
(0 until it.numRows).forEach { row -> it[row, 0] = get(row) }
|
||||
})
|
||||
}
|
||||
|
||||
public override fun buildMatrix(
|
||||
rows: Int,
|
||||
columns: Int,
|
||||
initializer: DoubleField.(i: Int, j: Int) -> Double,
|
||||
): EjmlDoubleMatrix<DMatrixRMaj> = EjmlDoubleMatrix(DMatrixRMaj(rows, columns).also {
|
||||
(0 until rows).forEach { row ->
|
||||
(0 until columns).forEach { col -> it[row, col] = elementAlgebra.initializer(row, col) }
|
||||
}
|
||||
})
|
||||
|
||||
public override fun buildVector(
|
||||
size: Int,
|
||||
initializer: DoubleField.(Int) -> Double,
|
||||
): EjmlDoubleVector<DMatrixRMaj> = EjmlDoubleVector(DMatrixRMaj(size, 1).also {
|
||||
(0 until it.numRows).forEach { row -> it[row, 0] = elementAlgebra.initializer(row) }
|
||||
})
|
||||
|
||||
private fun <T : DMatrix> T.wrapMatrix() = EjmlDoubleMatrix(this)
|
||||
private fun <T : DMatrixD1> T.wrapVector() = EjmlDoubleVector(this)
|
||||
|
||||
public override fun Matrix<Double>.unaryMinus(): Matrix<Double> = this * (-1.0)
|
||||
|
||||
public override fun Matrix<Double>.dot(other: Matrix<Double>): EjmlDoubleMatrix<DMatrixRMaj> {
|
||||
val out = DMatrixRMaj(1, 1)
|
||||
CommonOps_DDRM.mult(toEjml().origin, other.toEjml().origin, out)
|
||||
return out.wrapMatrix()
|
||||
}
|
||||
|
||||
public override fun Matrix<Double>.dot(vector: Point<Double>): EjmlDoubleVector<DMatrixRMaj> {
|
||||
val out = DMatrixRMaj(1, 1)
|
||||
CommonOps_DDRM.mult(toEjml().origin, vector.toEjml().origin, out)
|
||||
return out.wrapVector()
|
||||
}
|
||||
|
||||
public override operator fun Matrix<Double>.minus(other: Matrix<Double>): EjmlDoubleMatrix<DMatrixRMaj> {
|
||||
val out = DMatrixRMaj(1, 1)
|
||||
CommonOps_DDRM.subtract(toEjml().origin, other.toEjml().origin, out)
|
||||
return out.wrapMatrix()
|
||||
}
|
||||
|
||||
public override operator fun Matrix<Double>.times(value: Double): EjmlDoubleMatrix<DMatrixRMaj> {
|
||||
val res = this.toEjml().origin.copy()
|
||||
CommonOps_DDRM.scale(value, res)
|
||||
return res.wrapMatrix()
|
||||
}
|
||||
|
||||
public override fun Point<Double>.unaryMinus(): EjmlDoubleVector<DMatrixRMaj> {
|
||||
val out = toEjml().origin.copy()
|
||||
CommonOps_DDRM.changeSign(out)
|
||||
return out.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Matrix<Double>.plus(other: Matrix<Double>): EjmlDoubleMatrix<DMatrixRMaj> {
|
||||
val out = DMatrixRMaj(1, 1)
|
||||
CommonOps_DDRM.add(toEjml().origin, other.toEjml().origin, out)
|
||||
return out.wrapMatrix()
|
||||
}
|
||||
|
||||
public override fun Point<Double>.plus(other: Point<Double>): EjmlDoubleVector<DMatrixRMaj> {
|
||||
val out = DMatrixRMaj(1, 1)
|
||||
CommonOps_DDRM.add(toEjml().origin, other.toEjml().origin, out)
|
||||
return out.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Point<Double>.minus(other: Point<Double>): EjmlDoubleVector<DMatrixRMaj> {
|
||||
val out = DMatrixRMaj(1, 1)
|
||||
CommonOps_DDRM.subtract(toEjml().origin, other.toEjml().origin, out)
|
||||
return out.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Double.times(m: Matrix<Double>): EjmlDoubleMatrix<DMatrixRMaj> = m * this
|
||||
|
||||
public override fun Point<Double>.times(value: Double): EjmlDoubleVector<DMatrixRMaj> {
|
||||
val res = this.toEjml().origin.copy()
|
||||
CommonOps_DDRM.scale(value, res)
|
||||
return res.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Double.times(v: Point<Double>): EjmlDoubleVector<DMatrixRMaj> = v * this
|
||||
|
||||
@UnstableKMathAPI
|
||||
public override fun <F : StructureFeature> getFeature(structure: Matrix<Double>, type: KClass<out F>): F? {
|
||||
// Return the feature if it is intrinsic to the structure
|
||||
structure.getFeature(type)?.let { return it }
|
||||
|
||||
val origin = structure.toEjml().origin
|
||||
|
||||
return when (type) {
|
||||
InverseMatrixFeature::class -> object : InverseMatrixFeature<Double> {
|
||||
override val inverse: Matrix<Double> by lazy {
|
||||
val res = origin.copy()
|
||||
CommonOps_DDRM.invert(res)
|
||||
EjmlDoubleMatrix(res)
|
||||
}
|
||||
}
|
||||
|
||||
DeterminantFeature::class -> object : DeterminantFeature<Double> {
|
||||
override val determinant: Double by lazy { CommonOps_DDRM.det(DMatrixRMaj(origin)) }
|
||||
}
|
||||
|
||||
SingularValueDecompositionFeature::class -> object : SingularValueDecompositionFeature<Double> {
|
||||
private val svd by lazy {
|
||||
DecompositionFactory_DDRM.svd(origin.numRows, origin.numCols, true, true, false)
|
||||
.apply { decompose(origin.copy()) }
|
||||
}
|
||||
|
||||
override val u: Matrix<Double> by lazy { EjmlDoubleMatrix(svd.getU(null, false)) }
|
||||
override val s: Matrix<Double> by lazy { EjmlDoubleMatrix(svd.getW(null)) }
|
||||
override val v: Matrix<Double> by lazy { EjmlDoubleMatrix(svd.getV(null, false)) }
|
||||
override val singularValues: Point<Double> by lazy { DoubleBuffer(svd.singularValues) }
|
||||
}
|
||||
|
||||
QRDecompositionFeature::class -> object : QRDecompositionFeature<Double> {
|
||||
private val qr by lazy {
|
||||
DecompositionFactory_DDRM.qr().apply { decompose(origin.copy()) }
|
||||
}
|
||||
|
||||
override val q: Matrix<Double> by lazy {
|
||||
EjmlDoubleMatrix(qr.getQ(null, false)) + OrthogonalFeature
|
||||
}
|
||||
|
||||
override val r: Matrix<Double> by lazy { EjmlDoubleMatrix(qr.getR(null, false)) + UFeature }
|
||||
}
|
||||
|
||||
CholeskyDecompositionFeature::class -> object : CholeskyDecompositionFeature<Double> {
|
||||
override val l: Matrix<Double> by lazy {
|
||||
val cholesky =
|
||||
DecompositionFactory_DDRM.chol(structure.rowNum, true).apply { decompose(origin.copy()) }
|
||||
|
||||
EjmlDoubleMatrix(cholesky.getT(null)) + LFeature
|
||||
}
|
||||
}
|
||||
|
||||
LupDecompositionFeature::class -> object : LupDecompositionFeature<Double> {
|
||||
private val lup by lazy {
|
||||
DecompositionFactory_DDRM.lu(origin.numRows, origin.numCols).apply { decompose(origin.copy()) }
|
||||
}
|
||||
|
||||
override val l: Matrix<Double> by lazy {
|
||||
EjmlDoubleMatrix(lup.getLower(null)) + LFeature
|
||||
}
|
||||
|
||||
override val u: Matrix<Double> by lazy {
|
||||
EjmlDoubleMatrix(lup.getUpper(null)) + UFeature
|
||||
}
|
||||
|
||||
override val p: Matrix<Double> by lazy { EjmlDoubleMatrix(lup.getRowPivot(null)) }
|
||||
}
|
||||
|
||||
else -> null
|
||||
}?.let(type::cast)
|
||||
}
|
||||
|
||||
/**
|
||||
* Solves for *x* in the following equation: *x = [a] <sup>-1</sup> · [b]*.
|
||||
*
|
||||
* @param a the base matrix.
|
||||
* @param b n by p matrix.
|
||||
* @return the solution for 'x' that is n by p.
|
||||
*/
|
||||
public fun solve(a: Matrix<Double>, b: Matrix<Double>): EjmlDoubleMatrix<DMatrixRMaj> {
|
||||
val res = DMatrixRMaj(1, 1)
|
||||
CommonOps_DDRM.solve(DMatrixRMaj(a.toEjml().origin), DMatrixRMaj(b.toEjml().origin), res)
|
||||
return EjmlDoubleMatrix(res)
|
||||
}
|
||||
|
||||
/**
|
||||
* Solves for *x* in the following equation: *x = [a] <sup>-1</sup> · [b]*.
|
||||
*
|
||||
* @param a the base matrix.
|
||||
* @param b n by p vector.
|
||||
* @return the solution for 'x' that is n by p.
|
||||
*/
|
||||
public fun solve(a: Matrix<Double>, b: Point<Double>): EjmlDoubleVector<DMatrixRMaj> {
|
||||
val res = DMatrixRMaj(1, 1)
|
||||
CommonOps_DDRM.solve(DMatrixRMaj(a.toEjml().origin), DMatrixRMaj(b.toEjml().origin), res)
|
||||
return EjmlDoubleVector(res)
|
||||
}
|
||||
}
|
||||
|
||||
@ -5,7 +5,6 @@
|
||||
|
||||
package space.kscience.kmath.ejml
|
||||
|
||||
import org.ejml.data.DMatrix
|
||||
import org.ejml.data.Matrix
|
||||
import space.kscience.kmath.nd.Structure2D
|
||||
|
||||
@ -21,12 +20,3 @@ public abstract class EjmlMatrix<T, out M : Matrix>(public open val origin: M) :
|
||||
public override val rowNum: Int get() = origin.numRows
|
||||
public override val colNum: Int get() = origin.numCols
|
||||
}
|
||||
|
||||
/**
|
||||
* [EjmlMatrix] specialization for [Double].
|
||||
*
|
||||
* @author Iaroslav Postovalov
|
||||
*/
|
||||
public class EjmlDoubleMatrix<out M : DMatrix>(public override val origin: M) : EjmlMatrix<Double, M>(origin) {
|
||||
public override operator fun get(i: Int, j: Int): Double = origin[i, j]
|
||||
}
|
||||
|
||||
@ -5,7 +5,6 @@
|
||||
|
||||
package space.kscience.kmath.ejml
|
||||
|
||||
import org.ejml.data.DMatrixD1
|
||||
import org.ejml.data.Matrix
|
||||
import space.kscience.kmath.linear.Point
|
||||
|
||||
@ -14,12 +13,12 @@ import space.kscience.kmath.linear.Point
|
||||
*
|
||||
* @param T the type of elements contained in the buffer.
|
||||
* @param M the type of EJML matrix.
|
||||
* @property origin The underlying matrix.
|
||||
* @property origin The underlying matrix, must have only one row.
|
||||
* @author Iaroslav Postovalov
|
||||
*/
|
||||
public abstract class EjmlVector<out T, out M : Matrix>(public open val origin: M) : Point<T> {
|
||||
public override val size: Int
|
||||
get() = origin.numRows
|
||||
get() = origin.numCols
|
||||
|
||||
public override operator fun iterator(): Iterator<T> = object : Iterator<T> {
|
||||
private var cursor: Int = 0
|
||||
@ -34,12 +33,3 @@ public abstract class EjmlVector<out T, out M : Matrix>(public open val origin:
|
||||
|
||||
public override fun toString(): String = "EjmlVector(origin=$origin)"
|
||||
}
|
||||
|
||||
/**
|
||||
* [EjmlVector] specialization for [Double].
|
||||
*
|
||||
* @author Iaroslav Postovalov
|
||||
*/
|
||||
public class EjmlDoubleVector<out M : DMatrixD1>(public override val origin: M) : EjmlVector<Double, M>(origin) {
|
||||
public override operator fun get(index: Int): Double = origin[index]
|
||||
}
|
||||
|
||||
@ -0,0 +1,995 @@
|
||||
/*
|
||||
* Copyright 2018-2021 KMath contributors.
|
||||
* Use of this source code is governed by the Apache 2.0 license that can be found in the license/LICENSE.txt file.
|
||||
*/
|
||||
|
||||
/* This file is generated with buildSrc/src/main/kotlin/space/kscience/kmath/ejml/codegen/ejmlCodegen.kt */
|
||||
|
||||
package space.kscience.kmath.ejml
|
||||
|
||||
import org.ejml.data.*
|
||||
import org.ejml.dense.row.CommonOps_DDRM
|
||||
import org.ejml.dense.row.CommonOps_FDRM
|
||||
import org.ejml.dense.row.factory.DecompositionFactory_DDRM
|
||||
import org.ejml.dense.row.factory.DecompositionFactory_FDRM
|
||||
import org.ejml.sparse.FillReducing
|
||||
import org.ejml.sparse.csc.CommonOps_DSCC
|
||||
import org.ejml.sparse.csc.CommonOps_FSCC
|
||||
import org.ejml.sparse.csc.factory.DecompositionFactory_DSCC
|
||||
import org.ejml.sparse.csc.factory.DecompositionFactory_FSCC
|
||||
import org.ejml.sparse.csc.factory.LinearSolverFactory_DSCC
|
||||
import org.ejml.sparse.csc.factory.LinearSolverFactory_FSCC
|
||||
import space.kscience.kmath.linear.*
|
||||
import space.kscience.kmath.linear.Matrix
|
||||
import space.kscience.kmath.misc.UnstableKMathAPI
|
||||
import space.kscience.kmath.nd.StructureFeature
|
||||
import space.kscience.kmath.operations.DoubleField
|
||||
import space.kscience.kmath.operations.FloatField
|
||||
import space.kscience.kmath.operations.invoke
|
||||
import space.kscience.kmath.structures.DoubleBuffer
|
||||
import space.kscience.kmath.structures.FloatBuffer
|
||||
import kotlin.reflect.KClass
|
||||
import kotlin.reflect.cast
|
||||
|
||||
/**
|
||||
* [EjmlVector] specialization for [Double].
|
||||
*/
|
||||
public class EjmlDoubleVector<out M : DMatrix>(public override val origin: M) : EjmlVector<Double, M>(origin) {
|
||||
init {
|
||||
require(origin.numRows == 1) { "The origin matrix must have only one row to form a vector" }
|
||||
}
|
||||
|
||||
public override operator fun get(index: Int): Double = origin[0, index]
|
||||
}
|
||||
|
||||
/**
|
||||
* [EjmlVector] specialization for [Float].
|
||||
*/
|
||||
public class EjmlFloatVector<out M : FMatrix>(public override val origin: M) : EjmlVector<Float, M>(origin) {
|
||||
init {
|
||||
require(origin.numRows == 1) { "The origin matrix must have only one row to form a vector" }
|
||||
}
|
||||
|
||||
public override operator fun get(index: Int): Float = origin[0, index]
|
||||
}
|
||||
|
||||
/**
|
||||
* [EjmlMatrix] specialization for [Double].
|
||||
*/
|
||||
public class EjmlDoubleMatrix<out M : DMatrix>(public override val origin: M) : EjmlMatrix<Double, M>(origin) {
|
||||
public override operator fun get(i: Int, j: Int): Double = origin[i, j]
|
||||
}
|
||||
|
||||
/**
|
||||
* [EjmlMatrix] specialization for [Float].
|
||||
*/
|
||||
public class EjmlFloatMatrix<out M : FMatrix>(public override val origin: M) : EjmlMatrix<Float, M>(origin) {
|
||||
public override operator fun get(i: Int, j: Int): Float = origin[i, j]
|
||||
}
|
||||
|
||||
/**
|
||||
* [EjmlLinearSpace] implementation based on [CommonOps_DDRM], [DecompositionFactory_DDRM] operations and
|
||||
* [DMatrixRMaj] matrices.
|
||||
*/
|
||||
public object EjmlLinearSpaceDDRM : EjmlLinearSpace<Double, DoubleField, DMatrixRMaj>() {
|
||||
/**
|
||||
* The [DoubleField] reference.
|
||||
*/
|
||||
public override val elementAlgebra: DoubleField get() = DoubleField
|
||||
|
||||
@Suppress("UNCHECKED_CAST")
|
||||
public override fun Matrix<Double>.toEjml(): EjmlDoubleMatrix<DMatrixRMaj> = when {
|
||||
this is EjmlDoubleMatrix<*> && origin is DMatrixRMaj -> this as EjmlDoubleMatrix<DMatrixRMaj>
|
||||
else -> buildMatrix(rowNum, colNum) { i, j -> get(i, j) }
|
||||
}
|
||||
|
||||
@Suppress("UNCHECKED_CAST")
|
||||
public override fun Point<Double>.toEjml(): EjmlDoubleVector<DMatrixRMaj> = when {
|
||||
this is EjmlDoubleVector<*> && origin is DMatrixRMaj -> this as EjmlDoubleVector<DMatrixRMaj>
|
||||
else -> EjmlDoubleVector(DMatrixRMaj(size, 1).also {
|
||||
(0 until it.numRows).forEach { row -> it[row, 0] = get(row) }
|
||||
})
|
||||
}
|
||||
|
||||
public override fun buildMatrix(
|
||||
rows: Int,
|
||||
columns: Int,
|
||||
initializer: DoubleField.(i: Int, j: Int) -> Double,
|
||||
): EjmlDoubleMatrix<DMatrixRMaj> = DMatrixRMaj(rows, columns).also {
|
||||
(0 until rows).forEach { row ->
|
||||
(0 until columns).forEach { col -> it[row, col] = elementAlgebra.initializer(row, col) }
|
||||
}
|
||||
}.wrapMatrix()
|
||||
|
||||
public override fun buildVector(
|
||||
size: Int,
|
||||
initializer: DoubleField.(Int) -> Double,
|
||||
): EjmlDoubleVector<DMatrixRMaj> = EjmlDoubleVector(DMatrixRMaj(size, 1).also {
|
||||
(0 until it.numRows).forEach { row -> it[row, 0] = elementAlgebra.initializer(row) }
|
||||
})
|
||||
|
||||
private fun <T : DMatrix> T.wrapMatrix() = EjmlDoubleMatrix(this)
|
||||
private fun <T : DMatrix> T.wrapVector() = EjmlDoubleVector(this)
|
||||
|
||||
public override fun Matrix<Double>.unaryMinus(): Matrix<Double> = this * elementAlgebra { -one }
|
||||
|
||||
public override fun Matrix<Double>.dot(other: Matrix<Double>): EjmlDoubleMatrix<DMatrixRMaj> {
|
||||
val out = DMatrixRMaj(1, 1)
|
||||
CommonOps_DDRM.mult(toEjml().origin, other.toEjml().origin, out)
|
||||
return out.wrapMatrix()
|
||||
}
|
||||
|
||||
public override fun Matrix<Double>.dot(vector: Point<Double>): EjmlDoubleVector<DMatrixRMaj> {
|
||||
val out = DMatrixRMaj(1, 1)
|
||||
CommonOps_DDRM.mult(toEjml().origin, vector.toEjml().origin, out)
|
||||
return out.wrapVector()
|
||||
}
|
||||
|
||||
public override operator fun Matrix<Double>.minus(other: Matrix<Double>): EjmlDoubleMatrix<DMatrixRMaj> {
|
||||
val out = DMatrixRMaj(1, 1)
|
||||
|
||||
CommonOps_DDRM.add(
|
||||
elementAlgebra.one,
|
||||
toEjml().origin,
|
||||
elementAlgebra { -one },
|
||||
other.toEjml().origin,
|
||||
out,
|
||||
)
|
||||
|
||||
return out.wrapMatrix()
|
||||
}
|
||||
|
||||
public override operator fun Matrix<Double>.times(value: Double): EjmlDoubleMatrix<DMatrixRMaj> {
|
||||
val res = DMatrixRMaj(1, 1)
|
||||
CommonOps_DDRM.scale(value, toEjml().origin, res)
|
||||
return res.wrapMatrix()
|
||||
}
|
||||
|
||||
public override fun Point<Double>.unaryMinus(): EjmlDoubleVector<DMatrixRMaj> {
|
||||
val res = DMatrixRMaj(1, 1)
|
||||
CommonOps_DDRM.changeSign(toEjml().origin, res)
|
||||
return res.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Matrix<Double>.plus(other: Matrix<Double>): EjmlDoubleMatrix<DMatrixRMaj> {
|
||||
val out = DMatrixRMaj(1, 1)
|
||||
|
||||
CommonOps_DDRM.add(
|
||||
elementAlgebra.one,
|
||||
toEjml().origin,
|
||||
elementAlgebra.one,
|
||||
other.toEjml().origin,
|
||||
out,
|
||||
)
|
||||
|
||||
return out.wrapMatrix()
|
||||
}
|
||||
|
||||
public override fun Point<Double>.plus(other: Point<Double>): EjmlDoubleVector<DMatrixRMaj> {
|
||||
val out = DMatrixRMaj(1, 1)
|
||||
|
||||
CommonOps_DDRM.add(
|
||||
elementAlgebra.one,
|
||||
toEjml().origin,
|
||||
elementAlgebra.one,
|
||||
other.toEjml().origin,
|
||||
out,
|
||||
)
|
||||
|
||||
return out.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Point<Double>.minus(other: Point<Double>): EjmlDoubleVector<DMatrixRMaj> {
|
||||
val out = DMatrixRMaj(1, 1)
|
||||
|
||||
CommonOps_DDRM.add(
|
||||
elementAlgebra.one,
|
||||
toEjml().origin,
|
||||
elementAlgebra { -one },
|
||||
other.toEjml().origin,
|
||||
out,
|
||||
)
|
||||
|
||||
return out.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Double.times(m: Matrix<Double>): EjmlDoubleMatrix<DMatrixRMaj> = m * this
|
||||
|
||||
public override fun Point<Double>.times(value: Double): EjmlDoubleVector<DMatrixRMaj> {
|
||||
val res = DMatrixRMaj(1, 1)
|
||||
CommonOps_DDRM.scale(value, toEjml().origin, res)
|
||||
return res.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Double.times(v: Point<Double>): EjmlDoubleVector<DMatrixRMaj> = v * this
|
||||
|
||||
@UnstableKMathAPI
|
||||
public override fun <F : StructureFeature> getFeature(structure: Matrix<Double>, type: KClass<out F>): F? {
|
||||
structure.getFeature(type)?.let { return it }
|
||||
val origin = structure.toEjml().origin
|
||||
|
||||
return when (type) {
|
||||
InverseMatrixFeature::class -> object : InverseMatrixFeature<Double> {
|
||||
override val inverse: Matrix<Double> by lazy {
|
||||
val res = origin.copy()
|
||||
CommonOps_DDRM.invert(res)
|
||||
res.wrapMatrix()
|
||||
}
|
||||
}
|
||||
|
||||
DeterminantFeature::class -> object : DeterminantFeature<Double> {
|
||||
override val determinant: Double by lazy { CommonOps_DDRM.det(origin) }
|
||||
}
|
||||
|
||||
SingularValueDecompositionFeature::class -> object : SingularValueDecompositionFeature<Double> {
|
||||
private val svd by lazy {
|
||||
DecompositionFactory_DDRM.svd(origin.numRows, origin.numCols, true, true, false)
|
||||
.apply { decompose(origin.copy()) }
|
||||
}
|
||||
|
||||
override val u: Matrix<Double> by lazy { svd.getU(null, false).wrapMatrix() }
|
||||
override val s: Matrix<Double> by lazy { svd.getW(null).wrapMatrix() }
|
||||
override val v: Matrix<Double> by lazy { svd.getV(null, false).wrapMatrix() }
|
||||
override val singularValues: Point<Double> by lazy { DoubleBuffer(svd.singularValues) }
|
||||
}
|
||||
|
||||
QRDecompositionFeature::class -> object : QRDecompositionFeature<Double> {
|
||||
private val qr by lazy {
|
||||
DecompositionFactory_DDRM.qr().apply { decompose(origin.copy()) }
|
||||
}
|
||||
|
||||
override val q: Matrix<Double> by lazy {
|
||||
qr.getQ(null, false).wrapMatrix() + OrthogonalFeature
|
||||
}
|
||||
|
||||
override val r: Matrix<Double> by lazy { qr.getR(null, false).wrapMatrix() + UFeature }
|
||||
}
|
||||
|
||||
CholeskyDecompositionFeature::class -> object : CholeskyDecompositionFeature<Double> {
|
||||
override val l: Matrix<Double> by lazy {
|
||||
val cholesky =
|
||||
DecompositionFactory_DDRM.chol(structure.rowNum, true).apply { decompose(origin.copy()) }
|
||||
|
||||
cholesky.getT(null).wrapMatrix() + LFeature
|
||||
}
|
||||
}
|
||||
|
||||
LupDecompositionFeature::class -> object : LupDecompositionFeature<Double> {
|
||||
private val lup by lazy {
|
||||
DecompositionFactory_DDRM.lu(origin.numRows, origin.numCols).apply { decompose(origin.copy()) }
|
||||
}
|
||||
|
||||
override val l: Matrix<Double> by lazy {
|
||||
lup.getLower(null).wrapMatrix() + LFeature
|
||||
}
|
||||
|
||||
override val u: Matrix<Double> by lazy {
|
||||
lup.getUpper(null).wrapMatrix() + UFeature
|
||||
}
|
||||
|
||||
override val p: Matrix<Double> by lazy { lup.getRowPivot(null).wrapMatrix() }
|
||||
}
|
||||
|
||||
else -> null
|
||||
}?.let(type::cast)
|
||||
}
|
||||
|
||||
/**
|
||||
* Solves for *x* in the following equation: *x = [a] <sup>-1</sup> · [b]*.
|
||||
*
|
||||
* @param a the base matrix.
|
||||
* @param b n by p matrix.
|
||||
* @return the solution for *x* that is n by p.
|
||||
*/
|
||||
public fun solve(a: Matrix<Double>, b: Matrix<Double>): EjmlDoubleMatrix<DMatrixRMaj> {
|
||||
val res = DMatrixRMaj(1, 1)
|
||||
CommonOps_DDRM.solve(DMatrixRMaj(a.toEjml().origin), DMatrixRMaj(b.toEjml().origin), res)
|
||||
return res.wrapMatrix()
|
||||
}
|
||||
|
||||
/**
|
||||
* Solves for *x* in the following equation: *x = [a] <sup>-1</sup> · [b]*.
|
||||
*
|
||||
* @param a the base matrix.
|
||||
* @param b n by p vector.
|
||||
* @return the solution for *x* that is n by p.
|
||||
*/
|
||||
public fun solve(a: Matrix<Double>, b: Point<Double>): EjmlDoubleVector<DMatrixRMaj> {
|
||||
val res = DMatrixRMaj(1, 1)
|
||||
CommonOps_DDRM.solve(DMatrixRMaj(a.toEjml().origin), DMatrixRMaj(b.toEjml().origin), res)
|
||||
return EjmlDoubleVector(res)
|
||||
}
|
||||
}
|
||||
|
||||
/**
|
||||
* [EjmlLinearSpace] implementation based on [CommonOps_FDRM], [DecompositionFactory_FDRM] operations and
|
||||
* [FMatrixRMaj] matrices.
|
||||
*/
|
||||
public object EjmlLinearSpaceFDRM : EjmlLinearSpace<Float, FloatField, FMatrixRMaj>() {
|
||||
/**
|
||||
* The [FloatField] reference.
|
||||
*/
|
||||
public override val elementAlgebra: FloatField get() = FloatField
|
||||
|
||||
@Suppress("UNCHECKED_CAST")
|
||||
public override fun Matrix<Float>.toEjml(): EjmlFloatMatrix<FMatrixRMaj> = when {
|
||||
this is EjmlFloatMatrix<*> && origin is FMatrixRMaj -> this as EjmlFloatMatrix<FMatrixRMaj>
|
||||
else -> buildMatrix(rowNum, colNum) { i, j -> get(i, j) }
|
||||
}
|
||||
|
||||
@Suppress("UNCHECKED_CAST")
|
||||
public override fun Point<Float>.toEjml(): EjmlFloatVector<FMatrixRMaj> = when {
|
||||
this is EjmlFloatVector<*> && origin is FMatrixRMaj -> this as EjmlFloatVector<FMatrixRMaj>
|
||||
else -> EjmlFloatVector(FMatrixRMaj(size, 1).also {
|
||||
(0 until it.numRows).forEach { row -> it[row, 0] = get(row) }
|
||||
})
|
||||
}
|
||||
|
||||
public override fun buildMatrix(
|
||||
rows: Int,
|
||||
columns: Int,
|
||||
initializer: FloatField.(i: Int, j: Int) -> Float,
|
||||
): EjmlFloatMatrix<FMatrixRMaj> = FMatrixRMaj(rows, columns).also {
|
||||
(0 until rows).forEach { row ->
|
||||
(0 until columns).forEach { col -> it[row, col] = elementAlgebra.initializer(row, col) }
|
||||
}
|
||||
}.wrapMatrix()
|
||||
|
||||
public override fun buildVector(
|
||||
size: Int,
|
||||
initializer: FloatField.(Int) -> Float,
|
||||
): EjmlFloatVector<FMatrixRMaj> = EjmlFloatVector(FMatrixRMaj(size, 1).also {
|
||||
(0 until it.numRows).forEach { row -> it[row, 0] = elementAlgebra.initializer(row) }
|
||||
})
|
||||
|
||||
private fun <T : FMatrix> T.wrapMatrix() = EjmlFloatMatrix(this)
|
||||
private fun <T : FMatrix> T.wrapVector() = EjmlFloatVector(this)
|
||||
|
||||
public override fun Matrix<Float>.unaryMinus(): Matrix<Float> = this * elementAlgebra { -one }
|
||||
|
||||
public override fun Matrix<Float>.dot(other: Matrix<Float>): EjmlFloatMatrix<FMatrixRMaj> {
|
||||
val out = FMatrixRMaj(1, 1)
|
||||
CommonOps_FDRM.mult(toEjml().origin, other.toEjml().origin, out)
|
||||
return out.wrapMatrix()
|
||||
}
|
||||
|
||||
public override fun Matrix<Float>.dot(vector: Point<Float>): EjmlFloatVector<FMatrixRMaj> {
|
||||
val out = FMatrixRMaj(1, 1)
|
||||
CommonOps_FDRM.mult(toEjml().origin, vector.toEjml().origin, out)
|
||||
return out.wrapVector()
|
||||
}
|
||||
|
||||
public override operator fun Matrix<Float>.minus(other: Matrix<Float>): EjmlFloatMatrix<FMatrixRMaj> {
|
||||
val out = FMatrixRMaj(1, 1)
|
||||
|
||||
CommonOps_FDRM.add(
|
||||
elementAlgebra.one,
|
||||
toEjml().origin,
|
||||
elementAlgebra { -one },
|
||||
other.toEjml().origin,
|
||||
out,
|
||||
)
|
||||
|
||||
return out.wrapMatrix()
|
||||
}
|
||||
|
||||
public override operator fun Matrix<Float>.times(value: Float): EjmlFloatMatrix<FMatrixRMaj> {
|
||||
val res = FMatrixRMaj(1, 1)
|
||||
CommonOps_FDRM.scale(value, toEjml().origin, res)
|
||||
return res.wrapMatrix()
|
||||
}
|
||||
|
||||
public override fun Point<Float>.unaryMinus(): EjmlFloatVector<FMatrixRMaj> {
|
||||
val res = FMatrixRMaj(1, 1)
|
||||
CommonOps_FDRM.changeSign(toEjml().origin, res)
|
||||
return res.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Matrix<Float>.plus(other: Matrix<Float>): EjmlFloatMatrix<FMatrixRMaj> {
|
||||
val out = FMatrixRMaj(1, 1)
|
||||
|
||||
CommonOps_FDRM.add(
|
||||
elementAlgebra.one,
|
||||
toEjml().origin,
|
||||
elementAlgebra.one,
|
||||
other.toEjml().origin,
|
||||
out,
|
||||
)
|
||||
|
||||
return out.wrapMatrix()
|
||||
}
|
||||
|
||||
public override fun Point<Float>.plus(other: Point<Float>): EjmlFloatVector<FMatrixRMaj> {
|
||||
val out = FMatrixRMaj(1, 1)
|
||||
|
||||
CommonOps_FDRM.add(
|
||||
elementAlgebra.one,
|
||||
toEjml().origin,
|
||||
elementAlgebra.one,
|
||||
other.toEjml().origin,
|
||||
out,
|
||||
)
|
||||
|
||||
return out.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Point<Float>.minus(other: Point<Float>): EjmlFloatVector<FMatrixRMaj> {
|
||||
val out = FMatrixRMaj(1, 1)
|
||||
|
||||
CommonOps_FDRM.add(
|
||||
elementAlgebra.one,
|
||||
toEjml().origin,
|
||||
elementAlgebra { -one },
|
||||
other.toEjml().origin,
|
||||
out,
|
||||
)
|
||||
|
||||
return out.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Float.times(m: Matrix<Float>): EjmlFloatMatrix<FMatrixRMaj> = m * this
|
||||
|
||||
public override fun Point<Float>.times(value: Float): EjmlFloatVector<FMatrixRMaj> {
|
||||
val res = FMatrixRMaj(1, 1)
|
||||
CommonOps_FDRM.scale(value, toEjml().origin, res)
|
||||
return res.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Float.times(v: Point<Float>): EjmlFloatVector<FMatrixRMaj> = v * this
|
||||
|
||||
@UnstableKMathAPI
|
||||
public override fun <F : StructureFeature> getFeature(structure: Matrix<Float>, type: KClass<out F>): F? {
|
||||
structure.getFeature(type)?.let { return it }
|
||||
val origin = structure.toEjml().origin
|
||||
|
||||
return when (type) {
|
||||
InverseMatrixFeature::class -> object : InverseMatrixFeature<Float> {
|
||||
override val inverse: Matrix<Float> by lazy {
|
||||
val res = origin.copy()
|
||||
CommonOps_FDRM.invert(res)
|
||||
res.wrapMatrix()
|
||||
}
|
||||
}
|
||||
|
||||
DeterminantFeature::class -> object : DeterminantFeature<Float> {
|
||||
override val determinant: Float by lazy { CommonOps_FDRM.det(origin) }
|
||||
}
|
||||
|
||||
SingularValueDecompositionFeature::class -> object : SingularValueDecompositionFeature<Float> {
|
||||
private val svd by lazy {
|
||||
DecompositionFactory_FDRM.svd(origin.numRows, origin.numCols, true, true, false)
|
||||
.apply { decompose(origin.copy()) }
|
||||
}
|
||||
|
||||
override val u: Matrix<Float> by lazy { svd.getU(null, false).wrapMatrix() }
|
||||
override val s: Matrix<Float> by lazy { svd.getW(null).wrapMatrix() }
|
||||
override val v: Matrix<Float> by lazy { svd.getV(null, false).wrapMatrix() }
|
||||
override val singularValues: Point<Float> by lazy { FloatBuffer(svd.singularValues) }
|
||||
}
|
||||
|
||||
QRDecompositionFeature::class -> object : QRDecompositionFeature<Float> {
|
||||
private val qr by lazy {
|
||||
DecompositionFactory_FDRM.qr().apply { decompose(origin.copy()) }
|
||||
}
|
||||
|
||||
override val q: Matrix<Float> by lazy {
|
||||
qr.getQ(null, false).wrapMatrix() + OrthogonalFeature
|
||||
}
|
||||
|
||||
override val r: Matrix<Float> by lazy { qr.getR(null, false).wrapMatrix() + UFeature }
|
||||
}
|
||||
|
||||
CholeskyDecompositionFeature::class -> object : CholeskyDecompositionFeature<Float> {
|
||||
override val l: Matrix<Float> by lazy {
|
||||
val cholesky =
|
||||
DecompositionFactory_FDRM.chol(structure.rowNum, true).apply { decompose(origin.copy()) }
|
||||
|
||||
cholesky.getT(null).wrapMatrix() + LFeature
|
||||
}
|
||||
}
|
||||
|
||||
LupDecompositionFeature::class -> object : LupDecompositionFeature<Float> {
|
||||
private val lup by lazy {
|
||||
DecompositionFactory_FDRM.lu(origin.numRows, origin.numCols).apply { decompose(origin.copy()) }
|
||||
}
|
||||
|
||||
override val l: Matrix<Float> by lazy {
|
||||
lup.getLower(null).wrapMatrix() + LFeature
|
||||
}
|
||||
|
||||
override val u: Matrix<Float> by lazy {
|
||||
lup.getUpper(null).wrapMatrix() + UFeature
|
||||
}
|
||||
|
||||
override val p: Matrix<Float> by lazy { lup.getRowPivot(null).wrapMatrix() }
|
||||
}
|
||||
|
||||
else -> null
|
||||
}?.let(type::cast)
|
||||
}
|
||||
|
||||
/**
|
||||
* Solves for *x* in the following equation: *x = [a] <sup>-1</sup> · [b]*.
|
||||
*
|
||||
* @param a the base matrix.
|
||||
* @param b n by p matrix.
|
||||
* @return the solution for *x* that is n by p.
|
||||
*/
|
||||
public fun solve(a: Matrix<Float>, b: Matrix<Float>): EjmlFloatMatrix<FMatrixRMaj> {
|
||||
val res = FMatrixRMaj(1, 1)
|
||||
CommonOps_FDRM.solve(FMatrixRMaj(a.toEjml().origin), FMatrixRMaj(b.toEjml().origin), res)
|
||||
return res.wrapMatrix()
|
||||
}
|
||||
|
||||
/**
|
||||
* Solves for *x* in the following equation: *x = [a] <sup>-1</sup> · [b]*.
|
||||
*
|
||||
* @param a the base matrix.
|
||||
* @param b n by p vector.
|
||||
* @return the solution for *x* that is n by p.
|
||||
*/
|
||||
public fun solve(a: Matrix<Float>, b: Point<Float>): EjmlFloatVector<FMatrixRMaj> {
|
||||
val res = FMatrixRMaj(1, 1)
|
||||
CommonOps_FDRM.solve(FMatrixRMaj(a.toEjml().origin), FMatrixRMaj(b.toEjml().origin), res)
|
||||
return EjmlFloatVector(res)
|
||||
}
|
||||
}
|
||||
|
||||
/**
|
||||
* [EjmlLinearSpace] implementation based on [CommonOps_DSCC], [DecompositionFactory_DSCC] operations and
|
||||
* [DMatrixSparseCSC] matrices.
|
||||
*/
|
||||
public object EjmlLinearSpaceDSCC : EjmlLinearSpace<Double, DoubleField, DMatrixSparseCSC>() {
|
||||
/**
|
||||
* The [DoubleField] reference.
|
||||
*/
|
||||
public override val elementAlgebra: DoubleField get() = DoubleField
|
||||
|
||||
@Suppress("UNCHECKED_CAST")
|
||||
public override fun Matrix<Double>.toEjml(): EjmlDoubleMatrix<DMatrixSparseCSC> = when {
|
||||
this is EjmlDoubleMatrix<*> && origin is DMatrixSparseCSC -> this as EjmlDoubleMatrix<DMatrixSparseCSC>
|
||||
else -> buildMatrix(rowNum, colNum) { i, j -> get(i, j) }
|
||||
}
|
||||
|
||||
@Suppress("UNCHECKED_CAST")
|
||||
public override fun Point<Double>.toEjml(): EjmlDoubleVector<DMatrixSparseCSC> = when {
|
||||
this is EjmlDoubleVector<*> && origin is DMatrixSparseCSC -> this as EjmlDoubleVector<DMatrixSparseCSC>
|
||||
else -> EjmlDoubleVector(DMatrixSparseCSC(size, 1).also {
|
||||
(0 until it.numRows).forEach { row -> it[row, 0] = get(row) }
|
||||
})
|
||||
}
|
||||
|
||||
public override fun buildMatrix(
|
||||
rows: Int,
|
||||
columns: Int,
|
||||
initializer: DoubleField.(i: Int, j: Int) -> Double,
|
||||
): EjmlDoubleMatrix<DMatrixSparseCSC> = DMatrixSparseCSC(rows, columns).also {
|
||||
(0 until rows).forEach { row ->
|
||||
(0 until columns).forEach { col -> it[row, col] = elementAlgebra.initializer(row, col) }
|
||||
}
|
||||
}.wrapMatrix()
|
||||
|
||||
public override fun buildVector(
|
||||
size: Int,
|
||||
initializer: DoubleField.(Int) -> Double,
|
||||
): EjmlDoubleVector<DMatrixSparseCSC> = EjmlDoubleVector(DMatrixSparseCSC(size, 1).also {
|
||||
(0 until it.numRows).forEach { row -> it[row, 0] = elementAlgebra.initializer(row) }
|
||||
})
|
||||
|
||||
private fun <T : DMatrix> T.wrapMatrix() = EjmlDoubleMatrix(this)
|
||||
private fun <T : DMatrix> T.wrapVector() = EjmlDoubleVector(this)
|
||||
|
||||
public override fun Matrix<Double>.unaryMinus(): Matrix<Double> = this * elementAlgebra { -one }
|
||||
|
||||
public override fun Matrix<Double>.dot(other: Matrix<Double>): EjmlDoubleMatrix<DMatrixSparseCSC> {
|
||||
val out = DMatrixSparseCSC(1, 1)
|
||||
CommonOps_DSCC.mult(toEjml().origin, other.toEjml().origin, out)
|
||||
return out.wrapMatrix()
|
||||
}
|
||||
|
||||
public override fun Matrix<Double>.dot(vector: Point<Double>): EjmlDoubleVector<DMatrixSparseCSC> {
|
||||
val out = DMatrixSparseCSC(1, 1)
|
||||
CommonOps_DSCC.mult(toEjml().origin, vector.toEjml().origin, out)
|
||||
return out.wrapVector()
|
||||
}
|
||||
|
||||
public override operator fun Matrix<Double>.minus(other: Matrix<Double>): EjmlDoubleMatrix<DMatrixSparseCSC> {
|
||||
val out = DMatrixSparseCSC(1, 1)
|
||||
|
||||
CommonOps_DSCC.add(
|
||||
elementAlgebra.one,
|
||||
toEjml().origin,
|
||||
elementAlgebra { -one },
|
||||
other.toEjml().origin,
|
||||
out,
|
||||
null,
|
||||
null,
|
||||
)
|
||||
|
||||
return out.wrapMatrix()
|
||||
}
|
||||
|
||||
public override operator fun Matrix<Double>.times(value: Double): EjmlDoubleMatrix<DMatrixSparseCSC> {
|
||||
val res = DMatrixSparseCSC(1, 1)
|
||||
CommonOps_DSCC.scale(value, toEjml().origin, res)
|
||||
return res.wrapMatrix()
|
||||
}
|
||||
|
||||
public override fun Point<Double>.unaryMinus(): EjmlDoubleVector<DMatrixSparseCSC> {
|
||||
val res = DMatrixSparseCSC(1, 1)
|
||||
CommonOps_DSCC.changeSign(toEjml().origin, res)
|
||||
return res.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Matrix<Double>.plus(other: Matrix<Double>): EjmlDoubleMatrix<DMatrixSparseCSC> {
|
||||
val out = DMatrixSparseCSC(1, 1)
|
||||
|
||||
CommonOps_DSCC.add(
|
||||
elementAlgebra.one,
|
||||
toEjml().origin,
|
||||
elementAlgebra.one,
|
||||
other.toEjml().origin,
|
||||
out,
|
||||
null,
|
||||
null,
|
||||
)
|
||||
|
||||
return out.wrapMatrix()
|
||||
}
|
||||
|
||||
public override fun Point<Double>.plus(other: Point<Double>): EjmlDoubleVector<DMatrixSparseCSC> {
|
||||
val out = DMatrixSparseCSC(1, 1)
|
||||
|
||||
CommonOps_DSCC.add(
|
||||
elementAlgebra.one,
|
||||
toEjml().origin,
|
||||
elementAlgebra.one,
|
||||
other.toEjml().origin,
|
||||
out,
|
||||
null,
|
||||
null,
|
||||
)
|
||||
|
||||
return out.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Point<Double>.minus(other: Point<Double>): EjmlDoubleVector<DMatrixSparseCSC> {
|
||||
val out = DMatrixSparseCSC(1, 1)
|
||||
|
||||
CommonOps_DSCC.add(
|
||||
elementAlgebra.one,
|
||||
toEjml().origin,
|
||||
elementAlgebra { -one },
|
||||
other.toEjml().origin,
|
||||
out,
|
||||
null,
|
||||
null,
|
||||
)
|
||||
|
||||
return out.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Double.times(m: Matrix<Double>): EjmlDoubleMatrix<DMatrixSparseCSC> = m * this
|
||||
|
||||
public override fun Point<Double>.times(value: Double): EjmlDoubleVector<DMatrixSparseCSC> {
|
||||
val res = DMatrixSparseCSC(1, 1)
|
||||
CommonOps_DSCC.scale(value, toEjml().origin, res)
|
||||
return res.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Double.times(v: Point<Double>): EjmlDoubleVector<DMatrixSparseCSC> = v * this
|
||||
|
||||
@UnstableKMathAPI
|
||||
public override fun <F : StructureFeature> getFeature(structure: Matrix<Double>, type: KClass<out F>): F? {
|
||||
structure.getFeature(type)?.let { return it }
|
||||
val origin = structure.toEjml().origin
|
||||
|
||||
return when (type) {
|
||||
QRDecompositionFeature::class -> object : QRDecompositionFeature<Double> {
|
||||
private val qr by lazy {
|
||||
DecompositionFactory_DSCC.qr(FillReducing.NONE).apply { decompose(origin.copy()) }
|
||||
}
|
||||
|
||||
override val q: Matrix<Double> by lazy {
|
||||
qr.getQ(null, false).wrapMatrix() + OrthogonalFeature
|
||||
}
|
||||
|
||||
override val r: Matrix<Double> by lazy { qr.getR(null, false).wrapMatrix() + UFeature }
|
||||
}
|
||||
|
||||
CholeskyDecompositionFeature::class -> object : CholeskyDecompositionFeature<Double> {
|
||||
override val l: Matrix<Double> by lazy {
|
||||
val cholesky =
|
||||
DecompositionFactory_DSCC.cholesky().apply { decompose(origin.copy()) }
|
||||
|
||||
(cholesky.getT(null) as DMatrix).wrapMatrix() + LFeature
|
||||
}
|
||||
}
|
||||
|
||||
LUDecompositionFeature::class, DeterminantFeature::class, InverseMatrixFeature::class -> object :
|
||||
LUDecompositionFeature<Double>, DeterminantFeature<Double>, InverseMatrixFeature<Double> {
|
||||
private val lu by lazy {
|
||||
DecompositionFactory_DSCC.lu(FillReducing.NONE).apply { decompose(origin.copy()) }
|
||||
}
|
||||
|
||||
override val l: Matrix<Double> by lazy {
|
||||
lu.getLower(null).wrapMatrix() + LFeature
|
||||
}
|
||||
|
||||
override val u: Matrix<Double> by lazy {
|
||||
lu.getUpper(null).wrapMatrix() + UFeature
|
||||
}
|
||||
|
||||
override val inverse: Matrix<Double> by lazy {
|
||||
var a = origin
|
||||
val inverse = DMatrixRMaj(1, 1)
|
||||
val solver = LinearSolverFactory_DSCC.lu(FillReducing.NONE)
|
||||
if (solver.modifiesA()) a = a.copy()
|
||||
val i = CommonOps_DDRM.identity(a.numRows)
|
||||
solver.solve(i, inverse)
|
||||
inverse.wrapMatrix()
|
||||
}
|
||||
|
||||
override val determinant: Double by lazy { elementAlgebra.number(lu.computeDeterminant().real) }
|
||||
}
|
||||
|
||||
else -> null
|
||||
}?.let(type::cast)
|
||||
}
|
||||
|
||||
/**
|
||||
* Solves for *x* in the following equation: *x = [a] <sup>-1</sup> · [b]*.
|
||||
*
|
||||
* @param a the base matrix.
|
||||
* @param b n by p matrix.
|
||||
* @return the solution for *x* that is n by p.
|
||||
*/
|
||||
public fun solve(a: Matrix<Double>, b: Matrix<Double>): EjmlDoubleMatrix<DMatrixSparseCSC> {
|
||||
val res = DMatrixSparseCSC(1, 1)
|
||||
CommonOps_DSCC.solve(DMatrixSparseCSC(a.toEjml().origin), DMatrixSparseCSC(b.toEjml().origin), res)
|
||||
return res.wrapMatrix()
|
||||
}
|
||||
|
||||
/**
|
||||
* Solves for *x* in the following equation: *x = [a] <sup>-1</sup> · [b]*.
|
||||
*
|
||||
* @param a the base matrix.
|
||||
* @param b n by p vector.
|
||||
* @return the solution for *x* that is n by p.
|
||||
*/
|
||||
public fun solve(a: Matrix<Double>, b: Point<Double>): EjmlDoubleVector<DMatrixSparseCSC> {
|
||||
val res = DMatrixSparseCSC(1, 1)
|
||||
CommonOps_DSCC.solve(DMatrixSparseCSC(a.toEjml().origin), DMatrixSparseCSC(b.toEjml().origin), res)
|
||||
return EjmlDoubleVector(res)
|
||||
}
|
||||
}
|
||||
|
||||
/**
|
||||
* [EjmlLinearSpace] implementation based on [CommonOps_FSCC], [DecompositionFactory_FSCC] operations and
|
||||
* [FMatrixSparseCSC] matrices.
|
||||
*/
|
||||
public object EjmlLinearSpaceFSCC : EjmlLinearSpace<Float, FloatField, FMatrixSparseCSC>() {
|
||||
/**
|
||||
* The [FloatField] reference.
|
||||
*/
|
||||
public override val elementAlgebra: FloatField get() = FloatField
|
||||
|
||||
@Suppress("UNCHECKED_CAST")
|
||||
public override fun Matrix<Float>.toEjml(): EjmlFloatMatrix<FMatrixSparseCSC> = when {
|
||||
this is EjmlFloatMatrix<*> && origin is FMatrixSparseCSC -> this as EjmlFloatMatrix<FMatrixSparseCSC>
|
||||
else -> buildMatrix(rowNum, colNum) { i, j -> get(i, j) }
|
||||
}
|
||||
|
||||
@Suppress("UNCHECKED_CAST")
|
||||
public override fun Point<Float>.toEjml(): EjmlFloatVector<FMatrixSparseCSC> = when {
|
||||
this is EjmlFloatVector<*> && origin is FMatrixSparseCSC -> this as EjmlFloatVector<FMatrixSparseCSC>
|
||||
else -> EjmlFloatVector(FMatrixSparseCSC(size, 1).also {
|
||||
(0 until it.numRows).forEach { row -> it[row, 0] = get(row) }
|
||||
})
|
||||
}
|
||||
|
||||
public override fun buildMatrix(
|
||||
rows: Int,
|
||||
columns: Int,
|
||||
initializer: FloatField.(i: Int, j: Int) -> Float,
|
||||
): EjmlFloatMatrix<FMatrixSparseCSC> = FMatrixSparseCSC(rows, columns).also {
|
||||
(0 until rows).forEach { row ->
|
||||
(0 until columns).forEach { col -> it[row, col] = elementAlgebra.initializer(row, col) }
|
||||
}
|
||||
}.wrapMatrix()
|
||||
|
||||
public override fun buildVector(
|
||||
size: Int,
|
||||
initializer: FloatField.(Int) -> Float,
|
||||
): EjmlFloatVector<FMatrixSparseCSC> = EjmlFloatVector(FMatrixSparseCSC(size, 1).also {
|
||||
(0 until it.numRows).forEach { row -> it[row, 0] = elementAlgebra.initializer(row) }
|
||||
})
|
||||
|
||||
private fun <T : FMatrix> T.wrapMatrix() = EjmlFloatMatrix(this)
|
||||
private fun <T : FMatrix> T.wrapVector() = EjmlFloatVector(this)
|
||||
|
||||
public override fun Matrix<Float>.unaryMinus(): Matrix<Float> = this * elementAlgebra { -one }
|
||||
|
||||
public override fun Matrix<Float>.dot(other: Matrix<Float>): EjmlFloatMatrix<FMatrixSparseCSC> {
|
||||
val out = FMatrixSparseCSC(1, 1)
|
||||
CommonOps_FSCC.mult(toEjml().origin, other.toEjml().origin, out)
|
||||
return out.wrapMatrix()
|
||||
}
|
||||
|
||||
public override fun Matrix<Float>.dot(vector: Point<Float>): EjmlFloatVector<FMatrixSparseCSC> {
|
||||
val out = FMatrixSparseCSC(1, 1)
|
||||
CommonOps_FSCC.mult(toEjml().origin, vector.toEjml().origin, out)
|
||||
return out.wrapVector()
|
||||
}
|
||||
|
||||
public override operator fun Matrix<Float>.minus(other: Matrix<Float>): EjmlFloatMatrix<FMatrixSparseCSC> {
|
||||
val out = FMatrixSparseCSC(1, 1)
|
||||
|
||||
CommonOps_FSCC.add(
|
||||
elementAlgebra.one,
|
||||
toEjml().origin,
|
||||
elementAlgebra { -one },
|
||||
other.toEjml().origin,
|
||||
out,
|
||||
null,
|
||||
null,
|
||||
)
|
||||
|
||||
return out.wrapMatrix()
|
||||
}
|
||||
|
||||
public override operator fun Matrix<Float>.times(value: Float): EjmlFloatMatrix<FMatrixSparseCSC> {
|
||||
val res = FMatrixSparseCSC(1, 1)
|
||||
CommonOps_FSCC.scale(value, toEjml().origin, res)
|
||||
return res.wrapMatrix()
|
||||
}
|
||||
|
||||
public override fun Point<Float>.unaryMinus(): EjmlFloatVector<FMatrixSparseCSC> {
|
||||
val res = FMatrixSparseCSC(1, 1)
|
||||
CommonOps_FSCC.changeSign(toEjml().origin, res)
|
||||
return res.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Matrix<Float>.plus(other: Matrix<Float>): EjmlFloatMatrix<FMatrixSparseCSC> {
|
||||
val out = FMatrixSparseCSC(1, 1)
|
||||
|
||||
CommonOps_FSCC.add(
|
||||
elementAlgebra.one,
|
||||
toEjml().origin,
|
||||
elementAlgebra.one,
|
||||
other.toEjml().origin,
|
||||
out,
|
||||
null,
|
||||
null,
|
||||
)
|
||||
|
||||
return out.wrapMatrix()
|
||||
}
|
||||
|
||||
public override fun Point<Float>.plus(other: Point<Float>): EjmlFloatVector<FMatrixSparseCSC> {
|
||||
val out = FMatrixSparseCSC(1, 1)
|
||||
|
||||
CommonOps_FSCC.add(
|
||||
elementAlgebra.one,
|
||||
toEjml().origin,
|
||||
elementAlgebra.one,
|
||||
other.toEjml().origin,
|
||||
out,
|
||||
null,
|
||||
null,
|
||||
)
|
||||
|
||||
return out.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Point<Float>.minus(other: Point<Float>): EjmlFloatVector<FMatrixSparseCSC> {
|
||||
val out = FMatrixSparseCSC(1, 1)
|
||||
|
||||
CommonOps_FSCC.add(
|
||||
elementAlgebra.one,
|
||||
toEjml().origin,
|
||||
elementAlgebra { -one },
|
||||
other.toEjml().origin,
|
||||
out,
|
||||
null,
|
||||
null,
|
||||
)
|
||||
|
||||
return out.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Float.times(m: Matrix<Float>): EjmlFloatMatrix<FMatrixSparseCSC> = m * this
|
||||
|
||||
public override fun Point<Float>.times(value: Float): EjmlFloatVector<FMatrixSparseCSC> {
|
||||
val res = FMatrixSparseCSC(1, 1)
|
||||
CommonOps_FSCC.scale(value, toEjml().origin, res)
|
||||
return res.wrapVector()
|
||||
}
|
||||
|
||||
public override fun Float.times(v: Point<Float>): EjmlFloatVector<FMatrixSparseCSC> = v * this
|
||||
|
||||
@UnstableKMathAPI
|
||||
public override fun <F : StructureFeature> getFeature(structure: Matrix<Float>, type: KClass<out F>): F? {
|
||||
structure.getFeature(type)?.let { return it }
|
||||
val origin = structure.toEjml().origin
|
||||
|
||||
return when (type) {
|
||||
QRDecompositionFeature::class -> object : QRDecompositionFeature<Float> {
|
||||
private val qr by lazy {
|
||||
DecompositionFactory_FSCC.qr(FillReducing.NONE).apply { decompose(origin.copy()) }
|
||||
}
|
||||
|
||||
override val q: Matrix<Float> by lazy {
|
||||
qr.getQ(null, false).wrapMatrix() + OrthogonalFeature
|
||||
}
|
||||
|
||||
override val r: Matrix<Float> by lazy { qr.getR(null, false).wrapMatrix() + UFeature }
|
||||
}
|
||||
|
||||
CholeskyDecompositionFeature::class -> object : CholeskyDecompositionFeature<Float> {
|
||||
override val l: Matrix<Float> by lazy {
|
||||
val cholesky =
|
||||
DecompositionFactory_FSCC.cholesky().apply { decompose(origin.copy()) }
|
||||
|
||||
(cholesky.getT(null) as FMatrix).wrapMatrix() + LFeature
|
||||
}
|
||||
}
|
||||
|
||||
LUDecompositionFeature::class, DeterminantFeature::class, InverseMatrixFeature::class -> object :
|
||||
LUDecompositionFeature<Float>, DeterminantFeature<Float>, InverseMatrixFeature<Float> {
|
||||
private val lu by lazy {
|
||||
DecompositionFactory_FSCC.lu(FillReducing.NONE).apply { decompose(origin.copy()) }
|
||||
}
|
||||
|
||||
override val l: Matrix<Float> by lazy {
|
||||
lu.getLower(null).wrapMatrix() + LFeature
|
||||
}
|
||||
|
||||
override val u: Matrix<Float> by lazy {
|
||||
lu.getUpper(null).wrapMatrix() + UFeature
|
||||
}
|
||||
|
||||
override val inverse: Matrix<Float> by lazy {
|
||||
var a = origin
|
||||
val inverse = FMatrixRMaj(1, 1)
|
||||
val solver = LinearSolverFactory_FSCC.lu(FillReducing.NONE)
|
||||
if (solver.modifiesA()) a = a.copy()
|
||||
val i = CommonOps_FDRM.identity(a.numRows)
|
||||
solver.solve(i, inverse)
|
||||
inverse.wrapMatrix()
|
||||
}
|
||||
|
||||
override val determinant: Float by lazy { elementAlgebra.number(lu.computeDeterminant().real) }
|
||||
}
|
||||
|
||||
else -> null
|
||||
}?.let(type::cast)
|
||||
}
|
||||
|
||||
/**
|
||||
* Solves for *x* in the following equation: *x = [a] <sup>-1</sup> · [b]*.
|
||||
*
|
||||
* @param a the base matrix.
|
||||
* @param b n by p matrix.
|
||||
* @return the solution for *x* that is n by p.
|
||||
*/
|
||||
public fun solve(a: Matrix<Float>, b: Matrix<Float>): EjmlFloatMatrix<FMatrixSparseCSC> {
|
||||
val res = FMatrixSparseCSC(1, 1)
|
||||
CommonOps_FSCC.solve(FMatrixSparseCSC(a.toEjml().origin), FMatrixSparseCSC(b.toEjml().origin), res)
|
||||
return res.wrapMatrix()
|
||||
}
|
||||
|
||||
/**
|
||||
* Solves for *x* in the following equation: *x = [a] <sup>-1</sup> · [b]*.
|
||||
*
|
||||
* @param a the base matrix.
|
||||
* @param b n by p vector.
|
||||
* @return the solution for *x* that is n by p.
|
||||
*/
|
||||
public fun solve(a: Matrix<Float>, b: Point<Float>): EjmlFloatVector<FMatrixSparseCSC> {
|
||||
val res = FMatrixSparseCSC(1, 1)
|
||||
CommonOps_FSCC.solve(FMatrixSparseCSC(a.toEjml().origin), FMatrixSparseCSC(b.toEjml().origin), res)
|
||||
return EjmlFloatVector(res)
|
||||
}
|
||||
}
|
||||
|
||||
@ -18,7 +18,7 @@ internal class EjmlVectorTest {
|
||||
|
||||
private val randomMatrix: DMatrixRMaj
|
||||
get() {
|
||||
val d = DMatrixRMaj(random.nextInt(2, 100), 1)
|
||||
val d = DMatrixRMaj(1, random.nextInt(2, 100))
|
||||
RandomMatrices_DDRM.fillUniform(d, random.asJavaRandom())
|
||||
return d
|
||||
}
|
||||
@ -27,7 +27,7 @@ internal class EjmlVectorTest {
|
||||
fun size() {
|
||||
val m = randomMatrix
|
||||
val w = EjmlDoubleVector(m)
|
||||
assertEquals(m.numRows, w.size)
|
||||
assertEquals(m.numCols, w.size)
|
||||
}
|
||||
|
||||
@Test
|
||||
@ -43,7 +43,7 @@ internal class EjmlVectorTest {
|
||||
val w = EjmlDoubleVector(m)
|
||||
|
||||
assertEquals(
|
||||
m.iterator(true, 0, 0, m.numRows - 1, 0).asSequence().toList(),
|
||||
m.iterator(true, 0, 0, 0, m.numCols - 1).asSequence().toList(),
|
||||
w.iterator().asSequence().toList()
|
||||
)
|
||||
}
|
||||
|
||||
Loading…
Reference in New Issue
Block a user