From 0aa73cd48fd0aa517189987ecbb0f011114d9a3e Mon Sep 17 00:00:00 2001 From: Alexander Nozik Date: Sun, 14 Mar 2021 09:43:22 +0300 Subject: [PATCH] Vector space refactor (optimization) --- examples/build.gradle.kts | 8 ++ .../kscience/kmath/benchmarks/DotBenchmark.kt | 17 ++-- .../ExpressionsInterpretersBenchmark.kt | 4 +- ...Benchmark.kt => MatrixInverseBenchmark.kt} | 13 +-- .../kmath/benchmarks/ViktorLogBenchmark.kt | 4 +- .../space/kscience/kmath/ast/MstAlgebra.kt | 6 +- kmath-core/api/kmath-core.api | 44 +++++++--- .../kscience/kmath/linear/LinearSpace.kt | 20 ++--- .../kmath/linear/LinearSpaceOverNd.kt | 85 +++++++++++++++++++ .../kscience/kmath/nd/BufferNDAlgebra.kt | 8 +- .../space/kscience/kmath/nd/Structure1D.kt | 6 +- .../space/kscience/kmath/nd/Structure2D.kt | 25 +++--- .../kmath/structures/BufferAccessor2D.kt | 4 +- .../kmath/structures/bufferOperation.kt | 7 ++ .../kscience/kmath/linear/RealLUSolverTest.kt | 4 +- .../kscience/kmath/ejml/EjmlLinearSpace.kt | 8 +- .../kotlin/space/kscience/kmath/real/dot.kt | 2 +- 17 files changed, 183 insertions(+), 82 deletions(-) rename examples/src/benchmarks/kotlin/space/kscience/kmath/benchmarks/{LinearAlgebraBenchmark.kt => MatrixInverseBenchmark.kt} (75%) create mode 100644 kmath-core/src/commonMain/kotlin/space/kscience/kmath/linear/LinearSpaceOverNd.kt create mode 100644 kmath-core/src/commonMain/kotlin/space/kscience/kmath/structures/bufferOperation.kt diff --git a/examples/build.gradle.kts b/examples/build.gradle.kts index 77152dc0a..0301e4d67 100644 --- a/examples/build.gradle.kts +++ b/examples/build.gradle.kts @@ -92,6 +92,14 @@ benchmark { iterationTimeUnit = "ms" // time unity for iterationTime, default is seconds include("ExpressionsInterpretersBenchmark") } + + configurations.register("matrixInverse") { + warmups = 1 // number of warmup iterations + iterations = 3 // number of iterations + iterationTime = 500 // time in seconds per iteration + iterationTimeUnit = "ms" // time unity for iterationTime, default is seconds + include("MatrixInverseBenchmark") + } } kotlin.sourceSets.all { diff --git a/examples/src/benchmarks/kotlin/space/kscience/kmath/benchmarks/DotBenchmark.kt b/examples/src/benchmarks/kotlin/space/kscience/kmath/benchmarks/DotBenchmark.kt index 8fb0c284e..dbf373929 100644 --- a/examples/src/benchmarks/kotlin/space/kscience/kmath/benchmarks/DotBenchmark.kt +++ b/examples/src/benchmarks/kotlin/space/kscience/kmath/benchmarks/DotBenchmark.kt @@ -6,12 +6,9 @@ import kotlinx.benchmark.Scope import kotlinx.benchmark.State import space.kscience.kmath.commons.linear.CMLinearSpace import space.kscience.kmath.ejml.EjmlLinearSpace -import space.kscience.kmath.linear.BufferLinearSpace -import space.kscience.kmath.linear.Matrix -import space.kscience.kmath.linear.RealLinearSpace +import space.kscience.kmath.linear.LinearSpace +import space.kscience.kmath.linear.invoke import space.kscience.kmath.operations.RealField -import space.kscience.kmath.operations.invoke -import space.kscience.kmath.structures.Buffer import kotlin.random.Random @State(Scope.Benchmark) @@ -21,8 +18,8 @@ internal class DotBenchmark { const val dim = 1000 //creating invertible matrix - val matrix1 = Matrix.real(dim, dim) { i, j -> if (i <= j) random.nextDouble() else 0.0 } - val matrix2 = Matrix.real(dim, dim) { i, j -> if (i <= j) random.nextDouble() else 0.0 } + val matrix1 = LinearSpace.real.buildMatrix(dim, dim) { i, j -> if (i <= j) random.nextDouble() else 0.0 } + val matrix2 = LinearSpace.real.buildMatrix(dim, dim) { i, j -> if (i <= j) random.nextDouble() else 0.0 } val cmMatrix1 = CMLinearSpace { matrix1.toCM() } val cmMatrix2 = CMLinearSpace { matrix2.toCM() } @@ -33,7 +30,7 @@ internal class DotBenchmark { @Benchmark fun cmDot(blackhole: Blackhole) { - CMLinearSpace { + CMLinearSpace.run { blackhole.consume(cmMatrix1 dot cmMatrix2) } } @@ -54,14 +51,14 @@ internal class DotBenchmark { @Benchmark fun bufferedDot(blackhole: Blackhole) { - BufferLinearSpace(RealField, Buffer.Companion::real).invoke { + LinearSpace.auto(RealField).invoke { blackhole.consume(matrix1 dot matrix2) } } @Benchmark fun realDot(blackhole: Blackhole) { - RealLinearSpace { + LinearSpace.real { blackhole.consume(matrix1 dot matrix2) } } diff --git a/examples/src/benchmarks/kotlin/space/kscience/kmath/benchmarks/ExpressionsInterpretersBenchmark.kt b/examples/src/benchmarks/kotlin/space/kscience/kmath/benchmarks/ExpressionsInterpretersBenchmark.kt index 9c150f074..e5cfbf9f6 100644 --- a/examples/src/benchmarks/kotlin/space/kscience/kmath/benchmarks/ExpressionsInterpretersBenchmark.kt +++ b/examples/src/benchmarks/kotlin/space/kscience/kmath/benchmarks/ExpressionsInterpretersBenchmark.kt @@ -30,7 +30,7 @@ internal class ExpressionsInterpretersBenchmark { fun mstExpression(blackhole: Blackhole) { val expr = algebra.mstInField { val x = bindSymbol(x) - x * 2.0 + 2.0 / x - 16.0 + x * 2.0 + number(2.0) / x - 16.0 } invokeAndSum(expr, blackhole) @@ -40,7 +40,7 @@ internal class ExpressionsInterpretersBenchmark { fun asmExpression(blackhole: Blackhole) { val expr = algebra.mstInField { val x = bindSymbol(x) - x * 2.0 + 2.0 / x - 16.0 + x * 2.0 + number(2.0) / x - 16.0 }.compile() invokeAndSum(expr, blackhole) diff --git a/examples/src/benchmarks/kotlin/space/kscience/kmath/benchmarks/LinearAlgebraBenchmark.kt b/examples/src/benchmarks/kotlin/space/kscience/kmath/benchmarks/MatrixInverseBenchmark.kt similarity index 75% rename from examples/src/benchmarks/kotlin/space/kscience/kmath/benchmarks/LinearAlgebraBenchmark.kt rename to examples/src/benchmarks/kotlin/space/kscience/kmath/benchmarks/MatrixInverseBenchmark.kt index 85759d93c..7aa8ac975 100644 --- a/examples/src/benchmarks/kotlin/space/kscience/kmath/benchmarks/LinearAlgebraBenchmark.kt +++ b/examples/src/benchmarks/kotlin/space/kscience/kmath/benchmarks/MatrixInverseBenchmark.kt @@ -9,21 +9,22 @@ import space.kscience.kmath.commons.linear.inverse import space.kscience.kmath.ejml.EjmlLinearSpace import space.kscience.kmath.ejml.inverse import space.kscience.kmath.linear.LinearSpace -import space.kscience.kmath.linear.Matrix import space.kscience.kmath.linear.inverseWithLup -import space.kscience.kmath.linear.real +import space.kscience.kmath.linear.invoke import kotlin.random.Random @State(Scope.Benchmark) -internal class LinearAlgebraBenchmark { +internal class MatrixInverseBenchmark { companion object { val random = Random(1224) const val dim = 100 + private val space = LinearSpace.real + //creating invertible matrix - val u = Matrix.real(dim, dim) { i, j -> if (i <= j) random.nextDouble() else 0.0 } - val l = Matrix.real(dim, dim) { i, j -> if (i >= j) random.nextDouble() else 0.0 } - val matrix = l dot u + val u = space.buildMatrix(dim, dim) { i, j -> if (i <= j) random.nextDouble() else 0.0 } + val l = space.buildMatrix(dim, dim) { i, j -> if (i >= j) random.nextDouble() else 0.0 } + val matrix = space { l dot u } } @Benchmark diff --git a/examples/src/benchmarks/kotlin/space/kscience/kmath/benchmarks/ViktorLogBenchmark.kt b/examples/src/benchmarks/kotlin/space/kscience/kmath/benchmarks/ViktorLogBenchmark.kt index c48c86af9..0036b615c 100644 --- a/examples/src/benchmarks/kotlin/space/kscience/kmath/benchmarks/ViktorLogBenchmark.kt +++ b/examples/src/benchmarks/kotlin/space/kscience/kmath/benchmarks/ViktorLogBenchmark.kt @@ -15,7 +15,7 @@ import space.kscience.kmath.viktor.ViktorNDField internal class ViktorLogBenchmark { @Benchmark fun realFieldLog(blackhole: Blackhole) { - with(realField) { + with(realNdField) { val fortyTwo = produce { 42.0 } var res = one repeat(n) { res = ln(fortyTwo) } @@ -47,7 +47,7 @@ internal class ViktorLogBenchmark { // automatically build context most suited for given type. private val autoField = NDAlgebra.auto(RealField, dim, dim) - private val realField = NDAlgebra.real(dim, dim) + private val realNdField = NDAlgebra.real(dim, dim) private val viktorField = ViktorNDField(intArrayOf(dim, dim)) } } diff --git a/kmath-ast/src/commonMain/kotlin/space/kscience/kmath/ast/MstAlgebra.kt b/kmath-ast/src/commonMain/kotlin/space/kscience/kmath/ast/MstAlgebra.kt index 16fdfbeac..5ed39687b 100644 --- a/kmath-ast/src/commonMain/kotlin/space/kscience/kmath/ast/MstAlgebra.kt +++ b/kmath-ast/src/commonMain/kotlin/space/kscience/kmath/ast/MstAlgebra.kt @@ -54,7 +54,7 @@ public object MstRing : Ring, NumbersAddOperations, ScaleOperations, NumbersAddOperations, ScaleOperations< public override val one: MST.Numeric get() = MstRing.one - public override fun bindSymbol(value: String): MST.Symbolic = MstRing.bindSymbol(value) + public override fun bindSymbol(value: String): MST.Symbolic = MstAlgebra.bindSymbol(value) public override fun number(value: Number): MST.Numeric = MstRing.number(value) public override fun add(a: MST, b: MST): MST.Binary = MstRing.add(a, b) @@ -112,7 +112,7 @@ public object MstExtendedField : ExtendedField, NumericAlgebra { public override val zero: MST.Numeric get() = MstField.zero public override val one: MST.Numeric get() = MstField.one - public override fun bindSymbol(value: String): MST.Symbolic = MstField.bindSymbol(value) + public override fun bindSymbol(value: String): MST.Symbolic = MstAlgebra.bindSymbol(value) public override fun number(value: Number): MST.Numeric = MstRing.number(value) public override fun sin(arg: MST): MST.Unary = unaryOperationFunction(TrigonometricOperations.SIN_OPERATION)(arg) public override fun cos(arg: MST): MST.Unary = unaryOperationFunction(TrigonometricOperations.COS_OPERATION)(arg) diff --git a/kmath-core/api/kmath-core.api b/kmath-core/api/kmath-core.api index c1f8d5ec5..38da35d6c 100644 --- a/kmath-core/api/kmath-core.api +++ b/kmath-core/api/kmath-core.api @@ -496,7 +496,6 @@ public final class space/kscience/kmath/linear/LinearSpace$Companion { } public final class space/kscience/kmath/linear/LinearSpace$DefaultImpls { - public static fun buildVector (Lspace/kscience/kmath/linear/LinearSpace;ILkotlin/jvm/functions/Function2;)Lspace/kscience/kmath/structures/Buffer; public static fun dot (Lspace/kscience/kmath/linear/LinearSpace;Lspace/kscience/kmath/nd/Structure2D;Lspace/kscience/kmath/nd/Structure2D;)Lspace/kscience/kmath/nd/Structure2D; public static fun dot (Lspace/kscience/kmath/linear/LinearSpace;Lspace/kscience/kmath/nd/Structure2D;Lspace/kscience/kmath/structures/Buffer;)Lspace/kscience/kmath/structures/Buffer; public static fun minus (Lspace/kscience/kmath/linear/LinearSpace;Lspace/kscience/kmath/nd/Structure2D;Lspace/kscience/kmath/nd/Structure2D;)Lspace/kscience/kmath/nd/Structure2D; @@ -511,6 +510,29 @@ public final class space/kscience/kmath/linear/LinearSpace$DefaultImpls { public static fun unaryMinus (Lspace/kscience/kmath/linear/LinearSpace;Lspace/kscience/kmath/structures/Buffer;)Lspace/kscience/kmath/structures/Buffer; } +public final class space/kscience/kmath/linear/LinearSpaceKt { + public static final fun invoke (Lspace/kscience/kmath/linear/LinearSpace;Lkotlin/jvm/functions/Function1;)Ljava/lang/Object; +} + +public final class space/kscience/kmath/linear/LinearSpaceOverNd : space/kscience/kmath/linear/LinearSpace { + public fun (Lspace/kscience/kmath/operations/Ring;Lkotlin/jvm/functions/Function2;)V + public fun buildMatrix (IILkotlin/jvm/functions/Function3;)Lspace/kscience/kmath/nd/Structure2D; + public fun buildVector (ILkotlin/jvm/functions/Function2;)Lspace/kscience/kmath/structures/Buffer; + public fun dot (Lspace/kscience/kmath/nd/Structure2D;Lspace/kscience/kmath/nd/Structure2D;)Lspace/kscience/kmath/nd/Structure2D; + public fun dot (Lspace/kscience/kmath/nd/Structure2D;Lspace/kscience/kmath/structures/Buffer;)Lspace/kscience/kmath/structures/Buffer; + public fun getElementAlgebra ()Lspace/kscience/kmath/operations/Ring; + public fun minus (Lspace/kscience/kmath/nd/Structure2D;Lspace/kscience/kmath/nd/Structure2D;)Lspace/kscience/kmath/nd/Structure2D; + public fun minus (Lspace/kscience/kmath/structures/Buffer;Lspace/kscience/kmath/structures/Buffer;)Lspace/kscience/kmath/structures/Buffer; + public fun plus (Lspace/kscience/kmath/nd/Structure2D;Lspace/kscience/kmath/nd/Structure2D;)Lspace/kscience/kmath/nd/Structure2D; + public fun plus (Lspace/kscience/kmath/structures/Buffer;Lspace/kscience/kmath/structures/Buffer;)Lspace/kscience/kmath/structures/Buffer; + public fun times (Ljava/lang/Object;Lspace/kscience/kmath/nd/Structure2D;)Lspace/kscience/kmath/nd/Structure2D; + public fun times (Ljava/lang/Object;Lspace/kscience/kmath/structures/Buffer;)Lspace/kscience/kmath/structures/Buffer; + public fun times (Lspace/kscience/kmath/nd/Structure2D;Ljava/lang/Object;)Lspace/kscience/kmath/nd/Structure2D; + public fun times (Lspace/kscience/kmath/structures/Buffer;Ljava/lang/Object;)Lspace/kscience/kmath/structures/Buffer; + public fun unaryMinus (Lspace/kscience/kmath/nd/Structure2D;)Lspace/kscience/kmath/nd/Structure2D; + public fun unaryMinus (Lspace/kscience/kmath/structures/Buffer;)Lspace/kscience/kmath/structures/Buffer; +} + public final class space/kscience/kmath/linear/LupDecomposition : space/kscience/kmath/linear/DeterminantFeature, space/kscience/kmath/linear/LupDecompositionFeature { public fun (Lspace/kscience/kmath/linear/LinearSpace;Lspace/kscience/kmath/operations/Field;Lspace/kscience/kmath/nd/Structure2D;[IZ)V public final fun getContext ()Lspace/kscience/kmath/linear/LinearSpace; @@ -566,12 +588,12 @@ public final class space/kscience/kmath/linear/MatrixWrapper : space/kscience/km public fun get (II)Ljava/lang/Object; public fun get ([I)Ljava/lang/Object; public fun getColNum ()I - public fun getColumns ()Lspace/kscience/kmath/structures/Buffer; + public fun getColumns ()Ljava/util/List; public fun getDimension ()I public final fun getFeatures ()Ljava/util/Set; public final fun getOrigin ()Lspace/kscience/kmath/nd/Structure2D; public fun getRowNum ()I - public fun getRows ()Lspace/kscience/kmath/structures/Buffer; + public fun getRows ()Ljava/util/List; public fun getShape ()[I public fun hashCode ()I public fun toString ()Ljava/lang/String; @@ -622,11 +644,11 @@ public final class space/kscience/kmath/linear/VirtualMatrix : space/kscience/km public fun get (II)Ljava/lang/Object; public fun get ([I)Ljava/lang/Object; public fun getColNum ()I - public fun getColumns ()Lspace/kscience/kmath/structures/Buffer; + public fun getColumns ()Ljava/util/List; public fun getDimension ()I public final fun getGenerator ()Lkotlin/jvm/functions/Function2; public fun getRowNum ()I - public fun getRows ()Lspace/kscience/kmath/structures/Buffer; + public fun getRows ()Ljava/util/List; public fun getShape ()[I public fun hashCode ()I } @@ -678,11 +700,11 @@ public final class space/kscience/kmath/nd/BufferNDAlgebra$DefaultImpls { public final class space/kscience/kmath/nd/BufferNDAlgebraKt { public static final fun field (Lspace/kscience/kmath/nd/NDAlgebra$Companion;Lspace/kscience/kmath/operations/Field;Lkotlin/jvm/functions/Function2;[I)Lspace/kscience/kmath/nd/BufferedNDField; + public static final fun group (Lspace/kscience/kmath/nd/NDAlgebra$Companion;Lspace/kscience/kmath/operations/Group;Lkotlin/jvm/functions/Function2;[I)Lspace/kscience/kmath/nd/BufferedNDGroup; public static final fun ndField (Lspace/kscience/kmath/operations/Field;Lkotlin/jvm/functions/Function2;[ILkotlin/jvm/functions/Function1;)Ljava/lang/Object; + public static final fun ndGroup (Lspace/kscience/kmath/operations/Group;Lkotlin/jvm/functions/Function2;[ILkotlin/jvm/functions/Function1;)Ljava/lang/Object; public static final fun ndRing (Lspace/kscience/kmath/operations/Ring;Lkotlin/jvm/functions/Function2;[ILkotlin/jvm/functions/Function1;)Ljava/lang/Object; - public static final fun ndSpace (Lspace/kscience/kmath/operations/Group;Lkotlin/jvm/functions/Function2;[ILkotlin/jvm/functions/Function1;)Ljava/lang/Object; public static final fun ring (Lspace/kscience/kmath/nd/NDAlgebra$Companion;Lspace/kscience/kmath/operations/Ring;Lkotlin/jvm/functions/Function2;[I)Lspace/kscience/kmath/nd/BufferedNDRing; - public static final fun space (Lspace/kscience/kmath/nd/NDAlgebra$Companion;Lspace/kscience/kmath/operations/Group;Lkotlin/jvm/functions/Function2;[I)Lspace/kscience/kmath/nd/BufferedNDGroup; } public class space/kscience/kmath/nd/BufferedNDField : space/kscience/kmath/nd/BufferedNDRing, space/kscience/kmath/nd/NDField { @@ -1103,9 +1125,9 @@ public abstract interface class space/kscience/kmath/nd/Structure2D : space/ksci public abstract fun get (II)Ljava/lang/Object; public abstract fun get ([I)Ljava/lang/Object; public abstract fun getColNum ()I - public abstract fun getColumns ()Lspace/kscience/kmath/structures/Buffer; + public abstract fun getColumns ()Ljava/util/List; public abstract fun getRowNum ()I - public abstract fun getRows ()Lspace/kscience/kmath/structures/Buffer; + public abstract fun getRows ()Ljava/util/List; public abstract fun getShape ()[I } @@ -1115,9 +1137,9 @@ public final class space/kscience/kmath/nd/Structure2D$Companion { public final class space/kscience/kmath/nd/Structure2D$DefaultImpls { public static fun elements (Lspace/kscience/kmath/nd/Structure2D;)Lkotlin/sequences/Sequence; public static fun get (Lspace/kscience/kmath/nd/Structure2D;[I)Ljava/lang/Object; - public static fun getColumns (Lspace/kscience/kmath/nd/Structure2D;)Lspace/kscience/kmath/structures/Buffer; + public static fun getColumns (Lspace/kscience/kmath/nd/Structure2D;)Ljava/util/List; public static fun getDimension (Lspace/kscience/kmath/nd/Structure2D;)I - public static fun getRows (Lspace/kscience/kmath/nd/Structure2D;)Lspace/kscience/kmath/structures/Buffer; + public static fun getRows (Lspace/kscience/kmath/nd/Structure2D;)Ljava/util/List; public static fun getShape (Lspace/kscience/kmath/nd/Structure2D;)[I } diff --git a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/linear/LinearSpace.kt b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/linear/LinearSpace.kt index 4a9d2c9b3..b0331bdea 100644 --- a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/linear/LinearSpace.kt +++ b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/linear/LinearSpace.kt @@ -33,8 +33,7 @@ public interface LinearSpace> { /** * Produces a point compatible with matrix space (and possibly optimized for it). */ - public fun buildVector(size: Int, initializer: A.(Int) -> T): Vector = - buildMatrix(1, size) { _, j -> initializer(j) }.as1D() + public fun buildVector(size: Int, initializer: A.(Int) -> T): Vector public operator fun Matrix.unaryMinus(): Matrix = buildMatrix(rowNum, colNum) { i, j -> -get(i, j) @@ -154,7 +153,7 @@ public interface LinearSpace> { /** * Gets a feature from the matrix. This function may return some additional features to - * [space.kscience.kmath.nd.NDStructure.getFeature]. + * [group.kscience.kmath.nd.NDStructure.getFeature]. * * @param F the type of feature. * @param m the matrix. @@ -172,16 +171,7 @@ public interface LinearSpace> { public fun > buffered( algebra: A, bufferFactory: BufferFactory = Buffer.Companion::boxing, - ): LinearSpace = object : LinearSpace { - override val elementAlgebra: A = algebra - - override fun buildMatrix( - rows: Int, columns: Int, - initializer: A.(i: Int, j: Int) -> T, - ): Matrix = NDStructure.buffered(intArrayOf(rows, columns), bufferFactory) { (i, j) -> - algebra.initializer(i, j) - }.as2D() - } + ): LinearSpace = LinearSpaceOverNd(algebra,bufferFactory) public val real: LinearSpace = buffered(RealField, Buffer.Companion::real) @@ -193,9 +183,11 @@ public interface LinearSpace> { } } +public operator fun , R> LS.invoke(block: LS.() -> R): R = run(block) + /** * Gets a feature from the matrix. This function may return some additional features to - * [space.kscience.kmath.nd.NDStructure.getFeature]. + * [group.kscience.kmath.nd.NDStructure.getFeature]. * * @param T the type of items in the matrices. * @param M the type of operated matrices. diff --git a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/linear/LinearSpaceOverNd.kt b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/linear/LinearSpaceOverNd.kt new file mode 100644 index 000000000..6a82f37ec --- /dev/null +++ b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/linear/LinearSpaceOverNd.kt @@ -0,0 +1,85 @@ +package space.kscience.kmath.linear + +import space.kscience.kmath.nd.* +import space.kscience.kmath.operations.Ring +import space.kscience.kmath.operations.invoke +import space.kscience.kmath.structures.Buffer +import space.kscience.kmath.structures.BufferFactory +import space.kscience.kmath.structures.VirtualBuffer +import space.kscience.kmath.structures.indices + + +public class LinearSpaceOverNd>( + override val elementAlgebra: A, + private val bufferFactory: BufferFactory, +) : LinearSpace { + + private fun ndRing( + rows: Int, + cols: Int, + ): BufferedNDRing = NDAlgebra.ring(elementAlgebra, bufferFactory, rows, cols) + + override fun buildMatrix(rows: Int, columns: Int, initializer: A.(i: Int, j: Int) -> T): Matrix = + ndRing(rows, columns).produce { (i, j) -> elementAlgebra.initializer(i, j) }.as2D() + + override fun buildVector(size: Int, initializer: A.(Int) -> T): Vector = + bufferFactory(size) { elementAlgebra.initializer(it) } + + override fun Matrix.unaryMinus(): Matrix = ndRing(rowNum, colNum).run { + unwrap().map { -it }.as2D() + } + + override fun Matrix.plus(other: Matrix): Matrix = ndRing(rowNum, colNum).run { + require(shape.contentEquals(other.shape)) { "Shape mismatch on Matrix::plus. Expected $shape but found ${other.shape}" } + unwrap().plus(other.unwrap()).as2D() + } + + override fun Matrix.minus(other: Matrix): Matrix = ndRing(rowNum, colNum).run { + require(shape.contentEquals(other.shape)) { "Shape mismatch on Matrix::plus. Expected $shape but found ${other.shape}" } + unwrap().minus(other.unwrap()).as2D() + } + + private fun Buffer.linearize() = if (this is VirtualBuffer) { + buildVector(size) { get(it) } + } else { + this + } + + override fun Matrix.dot(other: Matrix): Matrix { + require(colNum == other.rowNum) { "Matrix dot operation dimension mismatch: ($rowNum, $colNum) x (${other.rowNum}, ${other.colNum})" } + return elementAlgebra { + val rows = this@dot.rows.map{it.linearize()} + val columns = other.columns.map { it.linearize() } + //TODO optimize buffers + buildMatrix(rowNum, other.colNum) { i, j -> + val r = rows[i] + val c = columns[j] + var res = zero + for (l in r.indices) { + res += r[l] * c[l] + } + res + } + } + } + + override fun Matrix.dot(vector: Vector): Vector { + require(colNum == vector.size) { "Matrix dot vector operation dimension mismatch: ($rowNum, $colNum) x (${vector.size})" } + return elementAlgebra { + val rows = this@dot.rows + //TODO optimize buffers + buildVector(rowNum) { i -> + val r = rows[i] + var res = zero + for (j in r.indices) { + res += r[j] * vector[j] + } + res + } + } + } + + override fun Matrix.times(value: T): Matrix = ndRing(rowNum, colNum).run { + unwrap().map { it * value }.as2D() + } +} \ No newline at end of file diff --git a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/BufferNDAlgebra.kt b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/BufferNDAlgebra.kt index 58d1ea7b1..bce3a0830 100644 --- a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/BufferNDAlgebra.kt +++ b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/BufferNDAlgebra.kt @@ -79,20 +79,20 @@ public open class BufferedNDField>( override fun scale(a: NDStructure, value: Double): NDStructure = a.map { it * value } } -// space factories -public fun > NDAlgebra.Companion.space( +// group factories +public fun > NDAlgebra.Companion.group( space: A, bufferFactory: BufferFactory, vararg shape: Int, ): BufferedNDGroup = BufferedNDGroup(shape, space, bufferFactory) -public inline fun , R> A.ndSpace( +public inline fun , R> A.ndGroup( noinline bufferFactory: BufferFactory, vararg shape: Int, action: BufferedNDGroup.() -> R, ): R { contract { callsInPlace(action, InvocationKind.EXACTLY_ONCE) } - return NDAlgebra.space(this, bufferFactory, *shape).run(action) + return NDAlgebra.group(this, bufferFactory, *shape).run(action) } //ring factories diff --git a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/Structure1D.kt b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/Structure1D.kt index c54bfeed9..1335a4933 100644 --- a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/Structure1D.kt +++ b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/Structure1D.kt @@ -45,14 +45,12 @@ private inline class Buffer1DWrapper(val buffer: Buffer) : Structure1D /** * Represent a [NDStructure] as [Structure1D]. Throw error in case of dimension mismatch */ -public fun NDStructure.as1D(): Structure1D = if (shape.size == 1) { +public fun NDStructure.as1D(): Structure1D = this as? Structure1D ?: if (shape.size == 1) { when (this) { - is Structure1DWrapper -> this is NDBuffer -> Buffer1DWrapper(this.buffer) else -> Structure1DWrapper(this) } -} else - error("Can't create 1d-structure from ${shape.size}d-structure") +} else error("Can't create 1d-structure from ${shape.size}d-structure") /** * Represent this buffer as 1D structure diff --git a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/Structure2D.kt b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/Structure2D.kt index 09c830646..e9f8234e5 100644 --- a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/Structure2D.kt +++ b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/Structure2D.kt @@ -26,14 +26,14 @@ public interface Structure2D : NDStructure { /** * The buffer of rows of this structure. It gets elements from the structure dynamically. */ - public val rows: Buffer> - get() = VirtualBuffer(rowNum) { i -> VirtualBuffer(colNum) { j -> get(i, j) } } + public val rows: List> + get() = List(rowNum) { i -> VirtualBuffer(colNum) { j -> get(i, j) } } /** * The buffer of columns of this structure. It gets elements from the structure dynamically. */ - public val columns: Buffer> - get() = VirtualBuffer(colNum) { j -> VirtualBuffer(rowNum) { i -> get(i, j) } } + public val columns: List> + get() = List(colNum) { j -> VirtualBuffer(rowNum) { i -> get(i, j) } } /** * Retrieves an element from the structure by two indices. @@ -75,20 +75,15 @@ private class Structure2DWrapper(val structure: NDStructure) : Structure2D override fun equals(other: Any?): Boolean = structure == other - override fun hashCode(): Int = structure.hashCode() + override fun hashCode(): Int = structure.hashCode() } /** * Represent a [NDStructure] as [Structure1D]. Throw error in case of dimension mismatch */ -public fun NDStructure.as2D(): Structure2D = if (shape.size == 2) - Structure2DWrapper(this) -else - error("Can't create 2d-structure from ${shape.size}d-structure") +public fun NDStructure.as2D(): Structure2D = this as? Structure2D ?: when (shape.size) { + 2 -> Structure2DWrapper(this) + else -> error("Can't create 2d-structure from ${shape.size}d-structure") +} -/** - * Alias for [Structure2D] with more familiar name. - * - * @param T the type of items in the matrix. - */ -public typealias Matrix = Structure2D +internal fun Structure2D.unwrap(): NDStructure = if (this is Structure2DWrapper) structure else this \ No newline at end of file diff --git a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/structures/BufferAccessor2D.kt b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/structures/BufferAccessor2D.kt index d9e37ebd8..0910b2034 100644 --- a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/structures/BufferAccessor2D.kt +++ b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/structures/BufferAccessor2D.kt @@ -13,10 +13,10 @@ internal class BufferAccessor2D( public val colNum: Int, val factory: MutableBufferFactory, ) { - public operator fun Buffer.get(i: Int, j: Int): T = get(i + colNum * j) + public operator fun Buffer.get(i: Int, j: Int): T = get(i*colNum + j) public operator fun MutableBuffer.set(i: Int, j: Int, value: T) { - set(i + colNum * j, value) + set(i*colNum + j, value) } public inline fun create(crossinline init: (i: Int, j: Int) -> T): MutableBuffer = diff --git a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/structures/bufferOperation.kt b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/structures/bufferOperation.kt new file mode 100644 index 000000000..459136631 --- /dev/null +++ b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/structures/bufferOperation.kt @@ -0,0 +1,7 @@ +package space.kscience.kmath.structures + + +public inline fun Buffer.map( + bufferFactory: BufferFactory = Buffer.Companion::auto, + crossinline block: (T) -> R, +): Buffer = bufferFactory(size) { block(get(it)) } \ No newline at end of file diff --git a/kmath-core/src/commonTest/kotlin/space/kscience/kmath/linear/RealLUSolverTest.kt b/kmath-core/src/commonTest/kotlin/space/kscience/kmath/linear/RealLUSolverTest.kt index 704f91998..fb90b5e11 100644 --- a/kmath-core/src/commonTest/kotlin/space/kscience/kmath/linear/RealLUSolverTest.kt +++ b/kmath-core/src/commonTest/kotlin/space/kscience/kmath/linear/RealLUSolverTest.kt @@ -19,13 +19,13 @@ class RealLUSolverTest { LinearSpace.real.run { val matrix = matrix(2,2)( 3.0, 1.0, - 1.0, 3.0 + 2.0, 3.0 ) val lup = lup(matrix) //Check determinant - assertEquals(8.0, lup.determinant) + assertEquals(7.0, lup.determinant) assertEquals(lup.p dot matrix, lup.l dot lup.u) } diff --git a/kmath-ejml/src/main/kotlin/space/kscience/kmath/ejml/EjmlLinearSpace.kt b/kmath-ejml/src/main/kotlin/space/kscience/kmath/ejml/EjmlLinearSpace.kt index 57ce977fd..2d3a928a6 100644 --- a/kmath-ejml/src/main/kotlin/space/kscience/kmath/ejml/EjmlLinearSpace.kt +++ b/kmath-ejml/src/main/kotlin/space/kscience/kmath/ejml/EjmlLinearSpace.kt @@ -5,14 +5,13 @@ import space.kscience.kmath.linear.* import space.kscience.kmath.misc.UnstableKMathAPI import space.kscience.kmath.nd.getFeature import space.kscience.kmath.operations.RealField -import space.kscience.kmath.operations.ScaleOperations /** * Represents context of basic operations operating with [EjmlMatrix]. * * @author Iaroslav Postovalov */ -public object EjmlLinearSpace : LinearSpace, ScaleOperations> { +public object EjmlLinearSpace : LinearSpace { override val elementAlgebra: RealField get() = RealField @@ -50,7 +49,7 @@ public object EjmlLinearSpace : LinearSpace, ScaleOperations< private fun SimpleMatrix.wrapMatrix() = EjmlMatrix(this) private fun SimpleMatrix.wrapVector() = EjmlVector(this) - override fun Matrix.unaryMinus(): Matrix = this * (-1) + override fun Matrix.unaryMinus(): Matrix = this * (-1.0) public override fun Matrix.dot(other: Matrix): EjmlMatrix = EjmlMatrix(toEjml().origin.mult(other.toEjml().origin)) @@ -61,9 +60,6 @@ public object EjmlLinearSpace : LinearSpace, ScaleOperations< public override operator fun Matrix.minus(other: Matrix): EjmlMatrix = (toEjml().origin - other.toEjml().origin).wrapMatrix() - public override fun scale(a: Matrix, value: Double): EjmlMatrix = - a.toEjml().origin.scale(value).wrapMatrix() - public override operator fun Matrix.times(value: Double): EjmlMatrix = toEjml().origin.scale(value).wrapMatrix() diff --git a/kmath-for-real/src/commonMain/kotlin/space/kscience/kmath/real/dot.kt b/kmath-for-real/src/commonMain/kotlin/space/kscience/kmath/real/dot.kt index 1a7c1213c..0bc5c111b 100644 --- a/kmath-for-real/src/commonMain/kotlin/space/kscience/kmath/real/dot.kt +++ b/kmath-for-real/src/commonMain/kotlin/space/kscience/kmath/real/dot.kt @@ -1,7 +1,7 @@ package space.kscience.kmath.real import space.kscience.kmath.linear.LinearSpace -import space.kscience.kmath.nd.Matrix +import space.kscience.kmath.linear.Matrix /**