Refactor/linear algebra #241

Merged
altavir merged 11 commits from refactor/linear-algebra into dev 2021-03-14 21:49:47 +03:00
17 changed files with 183 additions and 82 deletions
Showing only changes of commit 0aa73cd48f - Show all commits

View File

@ -92,6 +92,14 @@ benchmark {
iterationTimeUnit = "ms" // time unity for iterationTime, default is seconds iterationTimeUnit = "ms" // time unity for iterationTime, default is seconds
include("ExpressionsInterpretersBenchmark") 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 { kotlin.sourceSets.all {

View File

@ -6,12 +6,9 @@ import kotlinx.benchmark.Scope
import kotlinx.benchmark.State import kotlinx.benchmark.State
import space.kscience.kmath.commons.linear.CMLinearSpace import space.kscience.kmath.commons.linear.CMLinearSpace
import space.kscience.kmath.ejml.EjmlLinearSpace import space.kscience.kmath.ejml.EjmlLinearSpace
import space.kscience.kmath.linear.BufferLinearSpace import space.kscience.kmath.linear.LinearSpace
import space.kscience.kmath.linear.Matrix import space.kscience.kmath.linear.invoke
import space.kscience.kmath.linear.RealLinearSpace
import space.kscience.kmath.operations.RealField import space.kscience.kmath.operations.RealField
import space.kscience.kmath.operations.invoke
import space.kscience.kmath.structures.Buffer
import kotlin.random.Random import kotlin.random.Random
@State(Scope.Benchmark) @State(Scope.Benchmark)
@ -21,8 +18,8 @@ internal class DotBenchmark {
const val dim = 1000 const val dim = 1000
//creating invertible matrix //creating invertible matrix
val matrix1 = 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 = Matrix.real(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 cmMatrix1 = CMLinearSpace { matrix1.toCM() }
val cmMatrix2 = CMLinearSpace { matrix2.toCM() } val cmMatrix2 = CMLinearSpace { matrix2.toCM() }
@ -33,7 +30,7 @@ internal class DotBenchmark {
@Benchmark @Benchmark
fun cmDot(blackhole: Blackhole) { fun cmDot(blackhole: Blackhole) {
CMLinearSpace { CMLinearSpace.run {
blackhole.consume(cmMatrix1 dot cmMatrix2) blackhole.consume(cmMatrix1 dot cmMatrix2)
} }
} }
@ -54,14 +51,14 @@ internal class DotBenchmark {
@Benchmark @Benchmark
fun bufferedDot(blackhole: Blackhole) { fun bufferedDot(blackhole: Blackhole) {
BufferLinearSpace(RealField, Buffer.Companion::real).invoke { LinearSpace.auto(RealField).invoke {
blackhole.consume(matrix1 dot matrix2) blackhole.consume(matrix1 dot matrix2)
} }
} }
@Benchmark @Benchmark
fun realDot(blackhole: Blackhole) { fun realDot(blackhole: Blackhole) {
RealLinearSpace { LinearSpace.real {
blackhole.consume(matrix1 dot matrix2) blackhole.consume(matrix1 dot matrix2)
} }
} }

View File

@ -30,7 +30,7 @@ internal class ExpressionsInterpretersBenchmark {
fun mstExpression(blackhole: Blackhole) { fun mstExpression(blackhole: Blackhole) {
val expr = algebra.mstInField { val expr = algebra.mstInField {
val x = bindSymbol(x) val x = bindSymbol(x)
x * 2.0 + 2.0 / x - 16.0 x * 2.0 + number(2.0) / x - 16.0
} }
invokeAndSum(expr, blackhole) invokeAndSum(expr, blackhole)
@ -40,7 +40,7 @@ internal class ExpressionsInterpretersBenchmark {
fun asmExpression(blackhole: Blackhole) { fun asmExpression(blackhole: Blackhole) {
val expr = algebra.mstInField { val expr = algebra.mstInField {
val x = bindSymbol(x) val x = bindSymbol(x)
x * 2.0 + 2.0 / x - 16.0 x * 2.0 + number(2.0) / x - 16.0
}.compile() }.compile()
invokeAndSum(expr, blackhole) invokeAndSum(expr, blackhole)

View File

@ -9,21 +9,22 @@ import space.kscience.kmath.commons.linear.inverse
import space.kscience.kmath.ejml.EjmlLinearSpace import space.kscience.kmath.ejml.EjmlLinearSpace
import space.kscience.kmath.ejml.inverse import space.kscience.kmath.ejml.inverse
import space.kscience.kmath.linear.LinearSpace import space.kscience.kmath.linear.LinearSpace
import space.kscience.kmath.linear.Matrix
import space.kscience.kmath.linear.inverseWithLup import space.kscience.kmath.linear.inverseWithLup
import space.kscience.kmath.linear.real import space.kscience.kmath.linear.invoke
import kotlin.random.Random import kotlin.random.Random
@State(Scope.Benchmark) @State(Scope.Benchmark)
internal class LinearAlgebraBenchmark { internal class MatrixInverseBenchmark {
companion object { companion object {
val random = Random(1224) val random = Random(1224)
const val dim = 100 const val dim = 100
private val space = LinearSpace.real
//creating invertible matrix //creating invertible matrix
val u = Matrix.real(dim, dim) { i, j -> if (i <= j) random.nextDouble() else 0.0 } val u = space.buildMatrix(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 l = space.buildMatrix(dim, dim) { i, j -> if (i >= j) random.nextDouble() else 0.0 }
val matrix = l dot u val matrix = space { l dot u }
} }
@Benchmark @Benchmark

View File

@ -15,7 +15,7 @@ import space.kscience.kmath.viktor.ViktorNDField
internal class ViktorLogBenchmark { internal class ViktorLogBenchmark {
@Benchmark @Benchmark
fun realFieldLog(blackhole: Blackhole) { fun realFieldLog(blackhole: Blackhole) {
with(realField) { with(realNdField) {
val fortyTwo = produce { 42.0 } val fortyTwo = produce { 42.0 }
var res = one var res = one
repeat(n) { res = ln(fortyTwo) } repeat(n) { res = ln(fortyTwo) }
@ -47,7 +47,7 @@ internal class ViktorLogBenchmark {
// automatically build context most suited for given type. // automatically build context most suited for given type.
private val autoField = NDAlgebra.auto(RealField, dim, dim) 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)) private val viktorField = ViktorNDField(intArrayOf(dim, dim))
} }
} }

View File

@ -54,7 +54,7 @@ public object MstRing : Ring<MST>, NumbersAddOperations<MST>, ScaleOperations<MS
public override val one: MST.Numeric = number(1.0) public override val one: MST.Numeric = number(1.0)
public override fun number(value: Number): MST.Numeric = MstGroup.number(value) public override fun number(value: Number): MST.Numeric = MstGroup.number(value)
public override fun bindSymbol(value: String): MST.Symbolic = MstGroup.bindSymbol(value) public override fun bindSymbol(value: String): MST.Symbolic = MstAlgebra.bindSymbol(value)
public override fun add(a: MST, b: MST): MST.Binary = MstGroup.add(a, b) public override fun add(a: MST, b: MST): MST.Binary = MstGroup.add(a, b)
public override fun scale(a: MST, value: Double): MST.Binary = public override fun scale(a: MST, value: Double): MST.Binary =
@ -83,7 +83,7 @@ public object MstField : Field<MST>, NumbersAddOperations<MST>, ScaleOperations<
public override val one: MST.Numeric get() = MstRing.one 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 number(value: Number): MST.Numeric = MstRing.number(value)
public override fun add(a: MST, b: MST): MST.Binary = MstRing.add(a, b) public override fun add(a: MST, b: MST): MST.Binary = MstRing.add(a, b)
@ -112,7 +112,7 @@ public object MstExtendedField : ExtendedField<MST>, NumericAlgebra<MST> {
public override val zero: MST.Numeric get() = MstField.zero public override val zero: MST.Numeric get() = MstField.zero
public override val one: MST.Numeric get() = MstField.one 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 number(value: Number): MST.Numeric = MstRing.number(value)
public override fun sin(arg: MST): MST.Unary = unaryOperationFunction(TrigonometricOperations.SIN_OPERATION)(arg) 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) public override fun cos(arg: MST): MST.Unary = unaryOperationFunction(TrigonometricOperations.COS_OPERATION)(arg)

View File

@ -496,7 +496,6 @@ public final class space/kscience/kmath/linear/LinearSpace$Companion {
} }
public final class space/kscience/kmath/linear/LinearSpace$DefaultImpls { 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/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 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; 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 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 <init> (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 final class space/kscience/kmath/linear/LupDecomposition : space/kscience/kmath/linear/DeterminantFeature, space/kscience/kmath/linear/LupDecompositionFeature {
public fun <init> (Lspace/kscience/kmath/linear/LinearSpace;Lspace/kscience/kmath/operations/Field;Lspace/kscience/kmath/nd/Structure2D;[IZ)V public fun <init> (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; 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 (II)Ljava/lang/Object;
public fun get ([I)Ljava/lang/Object; public fun get ([I)Ljava/lang/Object;
public fun getColNum ()I public fun getColNum ()I
public fun getColumns ()Lspace/kscience/kmath/structures/Buffer; public fun getColumns ()Ljava/util/List;
public fun getDimension ()I public fun getDimension ()I
public final fun getFeatures ()Ljava/util/Set; public final fun getFeatures ()Ljava/util/Set;
public final fun getOrigin ()Lspace/kscience/kmath/nd/Structure2D; public final fun getOrigin ()Lspace/kscience/kmath/nd/Structure2D;
public fun getRowNum ()I public fun getRowNum ()I
public fun getRows ()Lspace/kscience/kmath/structures/Buffer; public fun getRows ()Ljava/util/List;
public fun getShape ()[I public fun getShape ()[I
public fun hashCode ()I public fun hashCode ()I
public fun toString ()Ljava/lang/String; 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 (II)Ljava/lang/Object;
public fun get ([I)Ljava/lang/Object; public fun get ([I)Ljava/lang/Object;
public fun getColNum ()I public fun getColNum ()I
public fun getColumns ()Lspace/kscience/kmath/structures/Buffer; public fun getColumns ()Ljava/util/List;
public fun getDimension ()I public fun getDimension ()I
public final fun getGenerator ()Lkotlin/jvm/functions/Function2; public final fun getGenerator ()Lkotlin/jvm/functions/Function2;
public fun getRowNum ()I public fun getRowNum ()I
public fun getRows ()Lspace/kscience/kmath/structures/Buffer; public fun getRows ()Ljava/util/List;
public fun getShape ()[I public fun getShape ()[I
public fun hashCode ()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 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 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 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 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 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 { 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 (II)Ljava/lang/Object;
public abstract fun get ([I)Ljava/lang/Object; public abstract fun get ([I)Ljava/lang/Object;
public abstract fun getColNum ()I 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 getRowNum ()I
public abstract fun getRows ()Lspace/kscience/kmath/structures/Buffer; public abstract fun getRows ()Ljava/util/List;
public abstract fun getShape ()[I 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 final class space/kscience/kmath/nd/Structure2D$DefaultImpls {
public static fun elements (Lspace/kscience/kmath/nd/Structure2D;)Lkotlin/sequences/Sequence; 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 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 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 public static fun getShape (Lspace/kscience/kmath/nd/Structure2D;)[I
} }

View File

@ -33,8 +33,7 @@ public interface LinearSpace<T : Any, out A : Ring<T>> {
/** /**
* Produces a point compatible with matrix space (and possibly optimized for it). * Produces a point compatible with matrix space (and possibly optimized for it).
*/ */
public fun buildVector(size: Int, initializer: A.(Int) -> T): Vector<T> = public fun buildVector(size: Int, initializer: A.(Int) -> T): Vector<T>
buildMatrix(1, size) { _, j -> initializer(j) }.as1D()
public operator fun Matrix<T>.unaryMinus(): Matrix<T> = buildMatrix(rowNum, colNum) { i, j -> public operator fun Matrix<T>.unaryMinus(): Matrix<T> = buildMatrix(rowNum, colNum) { i, j ->
-get(i, j) -get(i, j)
@ -154,7 +153,7 @@ public interface LinearSpace<T : Any, out A : Ring<T>> {
/** /**
* Gets a feature from the matrix. This function may return some additional features to * 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 F the type of feature.
* @param m the matrix. * @param m the matrix.
@ -172,16 +171,7 @@ public interface LinearSpace<T : Any, out A : Ring<T>> {
public fun <T : Any, A : Ring<T>> buffered( public fun <T : Any, A : Ring<T>> buffered(
algebra: A, algebra: A,
bufferFactory: BufferFactory<T> = Buffer.Companion::boxing, bufferFactory: BufferFactory<T> = Buffer.Companion::boxing,
): LinearSpace<T, A> = object : LinearSpace<T, A> { ): LinearSpace<T, A> = LinearSpaceOverNd(algebra,bufferFactory)
override val elementAlgebra: A = algebra
override fun buildMatrix(
rows: Int, columns: Int,
initializer: A.(i: Int, j: Int) -> T,
): Matrix<T> = NDStructure.buffered(intArrayOf(rows, columns), bufferFactory) { (i, j) ->
algebra.initializer(i, j)
}.as2D()
}
public val real: LinearSpace<Double, RealField> = buffered(RealField, Buffer.Companion::real) public val real: LinearSpace<Double, RealField> = buffered(RealField, Buffer.Companion::real)
@ -193,9 +183,11 @@ public interface LinearSpace<T : Any, out A : Ring<T>> {
} }
} }
public operator fun <LS : LinearSpace<*, *>, R> LS.invoke(block: LS.() -> R): R = run(block)
/** /**
* Gets a feature from the matrix. This function may return some additional features to * 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 T the type of items in the matrices.
* @param M the type of operated matrices. * @param M the type of operated matrices.

View File

@ -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<T : Any, A : Ring<T>>(
override val elementAlgebra: A,
private val bufferFactory: BufferFactory<T>,
) : LinearSpace<T, A> {
private fun ndRing(
rows: Int,
cols: Int,
): BufferedNDRing<T, A> = NDAlgebra.ring(elementAlgebra, bufferFactory, rows, cols)
override fun buildMatrix(rows: Int, columns: Int, initializer: A.(i: Int, j: Int) -> T): Matrix<T> =
ndRing(rows, columns).produce { (i, j) -> elementAlgebra.initializer(i, j) }.as2D()
override fun buildVector(size: Int, initializer: A.(Int) -> T): Vector<T> =
bufferFactory(size) { elementAlgebra.initializer(it) }
override fun Matrix<T>.unaryMinus(): Matrix<T> = ndRing(rowNum, colNum).run {
unwrap().map { -it }.as2D()
}
override fun Matrix<T>.plus(other: Matrix<T>): Matrix<T> = 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<T>.minus(other: Matrix<T>): Matrix<T> = 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<T>.linearize() = if (this is VirtualBuffer) {
buildVector(size) { get(it) }
} else {
this
}
override fun Matrix<T>.dot(other: Matrix<T>): Matrix<T> {
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<T>.dot(vector: Vector<T>): Vector<T> {
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<T>.times(value: T): Matrix<T> = ndRing(rowNum, colNum).run {
unwrap().map { it * value }.as2D()
}
}

View File

@ -79,20 +79,20 @@ public open class BufferedNDField<T, R : Field<T>>(
override fun scale(a: NDStructure<T>, value: Double): NDStructure<T> = a.map { it * value } override fun scale(a: NDStructure<T>, value: Double): NDStructure<T> = a.map { it * value }
} }
// space factories // group factories
public fun <T, A : Group<T>> NDAlgebra.Companion.space( public fun <T, A : Group<T>> NDAlgebra.Companion.group(
space: A, space: A,
bufferFactory: BufferFactory<T>, bufferFactory: BufferFactory<T>,
vararg shape: Int, vararg shape: Int,
): BufferedNDGroup<T, A> = BufferedNDGroup(shape, space, bufferFactory) ): BufferedNDGroup<T, A> = BufferedNDGroup(shape, space, bufferFactory)
public inline fun <T, A : Group<T>, R> A.ndSpace( public inline fun <T, A : Group<T>, R> A.ndGroup(
noinline bufferFactory: BufferFactory<T>, noinline bufferFactory: BufferFactory<T>,
vararg shape: Int, vararg shape: Int,
action: BufferedNDGroup<T, A>.() -> R, action: BufferedNDGroup<T, A>.() -> R,
): R { ): R {
contract { callsInPlace(action, InvocationKind.EXACTLY_ONCE) } contract { callsInPlace(action, InvocationKind.EXACTLY_ONCE) }
return NDAlgebra.space(this, bufferFactory, *shape).run(action) return NDAlgebra.group(this, bufferFactory, *shape).run(action)
} }
//ring factories //ring factories

View File

@ -45,14 +45,12 @@ private inline class Buffer1DWrapper<T>(val buffer: Buffer<T>) : Structure1D<T>
/** /**
* Represent a [NDStructure] as [Structure1D]. Throw error in case of dimension mismatch * Represent a [NDStructure] as [Structure1D]. Throw error in case of dimension mismatch
*/ */
public fun <T> NDStructure<T>.as1D(): Structure1D<T> = if (shape.size == 1) { public fun <T> NDStructure<T>.as1D(): Structure1D<T> = this as? Structure1D<T> ?: if (shape.size == 1) {
when (this) { when (this) {
is Structure1DWrapper -> this
is NDBuffer -> Buffer1DWrapper(this.buffer) is NDBuffer -> Buffer1DWrapper(this.buffer)
else -> Structure1DWrapper(this) else -> Structure1DWrapper(this)
} }
} else } else error("Can't create 1d-structure from ${shape.size}d-structure")
error("Can't create 1d-structure from ${shape.size}d-structure")
/** /**
* Represent this buffer as 1D structure * Represent this buffer as 1D structure

View File

@ -26,14 +26,14 @@ public interface Structure2D<T> : NDStructure<T> {
/** /**
* The buffer of rows of this structure. It gets elements from the structure dynamically. * The buffer of rows of this structure. It gets elements from the structure dynamically.
*/ */
public val rows: Buffer<Buffer<T>> public val rows: List<Buffer<T>>
get() = VirtualBuffer(rowNum) { i -> VirtualBuffer(colNum) { j -> get(i, j) } } get() = List(rowNum) { i -> VirtualBuffer(colNum) { j -> get(i, j) } }
/** /**
* The buffer of columns of this structure. It gets elements from the structure dynamically. * The buffer of columns of this structure. It gets elements from the structure dynamically.
*/ */
public val columns: Buffer<Buffer<T>> public val columns: List<Buffer<T>>
get() = VirtualBuffer(colNum) { j -> VirtualBuffer(rowNum) { i -> get(i, j) } } get() = List(colNum) { j -> VirtualBuffer(rowNum) { i -> get(i, j) } }
/** /**
* Retrieves an element from the structure by two indices. * Retrieves an element from the structure by two indices.
@ -81,14 +81,9 @@ private class Structure2DWrapper<T>(val structure: NDStructure<T>) : Structure2D
/** /**
* Represent a [NDStructure] as [Structure1D]. Throw error in case of dimension mismatch * Represent a [NDStructure] as [Structure1D]. Throw error in case of dimension mismatch
*/ */
public fun <T> NDStructure<T>.as2D(): Structure2D<T> = if (shape.size == 2) public fun <T> NDStructure<T>.as2D(): Structure2D<T> = this as? Structure2D<T> ?: when (shape.size) {
Structure2DWrapper(this) 2 -> Structure2DWrapper(this)
else else -> error("Can't create 2d-structure from ${shape.size}d-structure")
error("Can't create 2d-structure from ${shape.size}d-structure") }
/** internal fun <T> Structure2D<T>.unwrap(): NDStructure<T> = if (this is Structure2DWrapper) structure else this
* Alias for [Structure2D] with more familiar name.
*
* @param T the type of items in the matrix.
*/
public typealias Matrix<T> = Structure2D<T>

View File

@ -13,10 +13,10 @@ internal class BufferAccessor2D<T : Any>(
public val colNum: Int, public val colNum: Int,
val factory: MutableBufferFactory<T>, val factory: MutableBufferFactory<T>,
) { ) {
public operator fun Buffer<T>.get(i: Int, j: Int): T = get(i + colNum * j) public operator fun Buffer<T>.get(i: Int, j: Int): T = get(i*colNum + j)
public operator fun MutableBuffer<T>.set(i: Int, j: Int, value: T) { public operator fun MutableBuffer<T>.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<T> = public inline fun create(crossinline init: (i: Int, j: Int) -> T): MutableBuffer<T> =

View File

@ -0,0 +1,7 @@
package space.kscience.kmath.structures
public inline fun <T : Any, reified R : Any> Buffer<T>.map(
bufferFactory: BufferFactory<R> = Buffer.Companion::auto,
crossinline block: (T) -> R,
): Buffer<R> = bufferFactory(size) { block(get(it)) }

View File

@ -19,13 +19,13 @@ class RealLUSolverTest {
LinearSpace.real.run { LinearSpace.real.run {
val matrix = matrix(2,2)( val matrix = matrix(2,2)(
3.0, 1.0, 3.0, 1.0,
1.0, 3.0 2.0, 3.0
) )
val lup = lup(matrix) val lup = lup(matrix)
//Check determinant //Check determinant
assertEquals(8.0, lup.determinant) assertEquals(7.0, lup.determinant)
assertEquals(lup.p dot matrix, lup.l dot lup.u) assertEquals(lup.p dot matrix, lup.l dot lup.u)
} }

View File

@ -5,14 +5,13 @@ import space.kscience.kmath.linear.*
import space.kscience.kmath.misc.UnstableKMathAPI import space.kscience.kmath.misc.UnstableKMathAPI
import space.kscience.kmath.nd.getFeature import space.kscience.kmath.nd.getFeature
import space.kscience.kmath.operations.RealField import space.kscience.kmath.operations.RealField
import space.kscience.kmath.operations.ScaleOperations
/** /**
* Represents context of basic operations operating with [EjmlMatrix]. * Represents context of basic operations operating with [EjmlMatrix].
* *
* @author Iaroslav Postovalov * @author Iaroslav Postovalov
*/ */
public object EjmlLinearSpace : LinearSpace<Double, RealField>, ScaleOperations<Matrix<Double>> { public object EjmlLinearSpace : LinearSpace<Double, RealField> {
override val elementAlgebra: RealField get() = RealField override val elementAlgebra: RealField get() = RealField
@ -50,7 +49,7 @@ public object EjmlLinearSpace : LinearSpace<Double, RealField>, ScaleOperations<
private fun SimpleMatrix.wrapMatrix() = EjmlMatrix(this) private fun SimpleMatrix.wrapMatrix() = EjmlMatrix(this)
private fun SimpleMatrix.wrapVector() = EjmlVector(this) private fun SimpleMatrix.wrapVector() = EjmlVector(this)
override fun Matrix<Double>.unaryMinus(): Matrix<Double> = this * (-1) override fun Matrix<Double>.unaryMinus(): Matrix<Double> = this * (-1.0)
public override fun Matrix<Double>.dot(other: Matrix<Double>): EjmlMatrix = public override fun Matrix<Double>.dot(other: Matrix<Double>): EjmlMatrix =
EjmlMatrix(toEjml().origin.mult(other.toEjml().origin)) EjmlMatrix(toEjml().origin.mult(other.toEjml().origin))
@ -61,9 +60,6 @@ public object EjmlLinearSpace : LinearSpace<Double, RealField>, ScaleOperations<
public override operator fun Matrix<Double>.minus(other: Matrix<Double>): EjmlMatrix = public override operator fun Matrix<Double>.minus(other: Matrix<Double>): EjmlMatrix =
(toEjml().origin - other.toEjml().origin).wrapMatrix() (toEjml().origin - other.toEjml().origin).wrapMatrix()
public override fun scale(a: Matrix<Double>, value: Double): EjmlMatrix =
a.toEjml().origin.scale(value).wrapMatrix()
public override operator fun Matrix<Double>.times(value: Double): EjmlMatrix = public override operator fun Matrix<Double>.times(value: Double): EjmlMatrix =
toEjml().origin.scale(value).wrapMatrix() toEjml().origin.scale(value).wrapMatrix()

View File

@ -1,7 +1,7 @@
package space.kscience.kmath.real package space.kscience.kmath.real
import space.kscience.kmath.linear.LinearSpace import space.kscience.kmath.linear.LinearSpace
import space.kscience.kmath.nd.Matrix import space.kscience.kmath.linear.Matrix
/** /**