From 5846f42141a210bef8ba873030f9566ebc14c5e9 Mon Sep 17 00:00:00 2001 From: Alexander Nozik Date: Fri, 15 Jul 2022 15:21:49 +0300 Subject: [PATCH] Grand derivative refactoring. Phase 1 --- .../kscience/kmath/expressions/DSCompiler.kt | 90 +- .../kmath/expressions/DerivativeStructure.kt | 191 ++-- .../DerivativeStructureExpression.kt | 327 +++--- .../kscience/kmath/linear/LinearSpace.kt | 2 +- .../space/kscience/kmath/misc/annotations.kt | 2 +- .../space/kscience/kmath/nd/BufferND.kt | 2 +- .../space/kscience/kmath/nd/StructureND.kt | 4 +- .../kmath/operations/DoubleBufferOps.kt | 4 +- .../kmath/operations/bufferOperation.kt | 26 +- .../space/kscience/kmath/structures/Buffer.kt | 8 +- .../space/kscience/kmath/ejml/_generated.kt | 1003 ----------------- .../histogram/UniformHistogramGroupND.kt | 6 +- .../kmath/multik/MultikDoubleAlgebra.kt | 7 + .../space/kscience/kmath/stat/Sampler.kt | 2 +- 14 files changed, 311 insertions(+), 1363 deletions(-) delete mode 100644 kmath-ejml/src/main/kotlin/space/kscience/kmath/ejml/_generated.kt diff --git a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/expressions/DSCompiler.kt b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/expressions/DSCompiler.kt index bb88ce52c..e0050cf03 100644 --- a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/expressions/DSCompiler.kt +++ b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/expressions/DSCompiler.kt @@ -5,11 +5,11 @@ package space.kscience.kmath.expressions + import space.kscience.kmath.operations.* import space.kscience.kmath.structures.Buffer import space.kscience.kmath.structures.MutableBuffer import space.kscience.kmath.structures.MutableBufferFactory -import kotlin.math.max import kotlin.math.min internal fun MutableBuffer.fill(element: T, fromIndex: Int = 0, toIndex: Int = size) { @@ -54,7 +54,7 @@ internal fun MutableBuffer.fill(element: T, fromIndex: Int = 0, toIndex: * @property order Derivation order. * @see DerivativeStructure */ -internal class DSCompiler> internal constructor( +class DSCompiler> internal constructor( val algebra: A, val bufferFactory: MutableBufferFactory, val freeParameters: Int, @@ -120,8 +120,7 @@ internal class DSCompiler> internal constructor( * This number includes the single 0 order derivative element, which is * guaranteed to be stored in the first element of the array. */ - val size: Int - get() = sizes[freeParameters][order] + val size: Int get() = sizes[freeParameters][order] /** * Get the index of a partial derivative in the array. @@ -178,7 +177,7 @@ internal fun DSCompiler.ln( operand: Buffer, operandOffset: Int, result: MutableBuffer, - resultOffset: Int + resultOffset: Int, ) where A : Field, A : ExponentialOperations = algebra { // create the function value and derivatives val function = bufferFactory(1 + order) { zero } @@ -211,7 +210,7 @@ internal fun DSCompiler.pow( operandOffset: Int, n: Int, result: MutableBuffer, - resultOffset: Int + resultOffset: Int, ) where A : Field, A : PowerOperations = algebra { if (n == 0) { // special case, x^0 = 1 for all x @@ -267,7 +266,7 @@ internal fun DSCompiler.exp( operand: Buffer, operandOffset: Int, result: MutableBuffer, - resultOffset: Int + resultOffset: Int, ) where A : Ring, A : ScaleOperations, A : ExponentialOperations = algebra { // create the function value and derivatives val function = bufferFactory(1 + order) { zero } @@ -290,7 +289,7 @@ internal fun DSCompiler.sqrt( operand: Buffer, operandOffset: Int, result: MutableBuffer, - resultOffset: Int + resultOffset: Int, ) where A : Field, A : PowerOperations = algebra { // create the function value and derivatives // [x^(1/n), (1/n)x^((1/n)-1), (1-n)/n^2x^((1/n)-2), ... ] @@ -351,7 +350,7 @@ internal fun DSCompiler.pow( operandOffset: Int, p: Double, result: MutableBuffer, - resultOffset: Int + resultOffset: Int, ) where A : Ring, A : NumericAlgebra, A : PowerOperations, A : ScaleOperations = algebra { // create the function value and derivatives // [x^p, px^(p-1), p(p-1)x^(p-2), ... ] @@ -387,7 +386,7 @@ internal fun DSCompiler.tan( operand: Buffer, operandOffset: Int, result: MutableBuffer, - resultOffset: Int + resultOffset: Int, ) where A : Ring, A : TrigonometricOperations, A : ScaleOperations = algebra { // create the function value and derivatives val function = bufferFactory(1 + order) { zero } @@ -469,7 +468,7 @@ internal fun DSCompiler.sin( operand: Buffer, operandOffset: Int, result: MutableBuffer, - resultOffset: Int + resultOffset: Int, ) where A : Ring, A : ScaleOperations, A : TrigonometricOperations = algebra { // create the function value and derivatives val function = bufferFactory(1 + order) { zero } @@ -497,7 +496,7 @@ internal fun DSCompiler.acos( operand: Buffer, operandOffset: Int, result: MutableBuffer, - resultOffset: Int + resultOffset: Int, ) where A : Field, A : TrigonometricOperations, A : PowerOperations = algebra { // create the function value and derivatives val function = bufferFactory(1 + order) { zero } @@ -559,7 +558,7 @@ internal fun DSCompiler.asin( operand: Buffer, operandOffset: Int, result: MutableBuffer, - resultOffset: Int + resultOffset: Int, ) where A : Field, A : TrigonometricOperations, A : PowerOperations = algebra { // create the function value and derivatives val function = bufferFactory(1 + order) { zero } @@ -618,7 +617,7 @@ internal fun DSCompiler.atan( operand: Buffer, operandOffset: Int, result: MutableBuffer, - resultOffset: Int + resultOffset: Int, ) where A : Field, A : TrigonometricOperations = algebra { // create the function value and derivatives val function = bufferFactory(1 + order) { zero } @@ -678,7 +677,7 @@ internal fun DSCompiler.cosh( operand: Buffer, operandOffset: Int, result: MutableBuffer, - resultOffset: Int + resultOffset: Int, ) where A : Ring, A : ScaleOperations, A : ExponentialOperations = algebra { // create the function value and derivatives val function = bufferFactory(1 + order) { zero } @@ -708,7 +707,7 @@ internal fun DSCompiler.tanh( operand: Buffer, operandOffset: Int, result: MutableBuffer, - resultOffset: Int + resultOffset: Int, ) where A : Field, A : ExponentialOperations = algebra { // create the function value and derivatives val function = bufferFactory(1 + order) { zero } @@ -765,7 +764,7 @@ internal fun DSCompiler.acosh( operand: Buffer, operandOffset: Int, result: MutableBuffer, - resultOffset: Int + resultOffset: Int, ) where A : Field, A : ExponentialOperations, A : PowerOperations = algebra { // create the function value and derivatives val function = bufferFactory(1 + order) { zero } @@ -857,7 +856,7 @@ internal fun DSCompiler.sinh( operand: Buffer, operandOffset: Int, result: MutableBuffer, - resultOffset: Int + resultOffset: Int, ) where A : Field, A : ExponentialOperations = algebra { // create the function value and derivatives val function = bufferFactory(1 + order) { zero } @@ -964,7 +963,7 @@ internal fun DSCompiler.asinh( operand: Buffer, operandOffset: Int, result: MutableBuffer, - resultOffset: Int + resultOffset: Int, ) where A : Field, A : ExponentialOperations, A : PowerOperations = algebra { // create the function value and derivatives val function = bufferFactory(1 + order) { zero } @@ -1109,59 +1108,6 @@ internal fun DSCompiler.atanh( compose(operand, operandOffset, function, result, resultOffset) } -/** - * Get the compiler for number of free parameters and order. - * - * @param parameters number of free parameters. - * @param order derivation order. - * @return cached rules set. - */ -internal fun > getCompiler( - algebra: A, - bufferFactory: MutableBufferFactory, - parameters: Int, - order: Int -): DSCompiler { - // get the cached compilers - val cache: Array?>>? = null - - // we need to create more compilers - val maxParameters: Int = max(parameters, cache?.size ?: 0) - val maxOrder: Int = max(order, if (cache == null) 0 else cache[0].size) - val newCache: Array?>> = Array(maxParameters + 1) { arrayOfNulls(maxOrder + 1) } - - if (cache != null) { - // preserve the already created compilers - for (i in cache.indices) { - cache[i].copyInto(newCache[i], endIndex = cache[i].size) - } - } - - // create the array in increasing diagonal order - - // create the array in increasing diagonal order - for (diag in 0..parameters + order) { - for (o in max(0, diag - parameters)..min(order, diag)) { - val p: Int = diag - o - if (newCache[p][o] == null) { - val valueCompiler: DSCompiler? = if (p == 0) null else newCache[p - 1][o]!! - val derivativeCompiler: DSCompiler? = if (o == 0) null else newCache[p][o - 1]!! - - newCache[p][o] = DSCompiler( - algebra, - bufferFactory, - p, - o, - valueCompiler, - derivativeCompiler, - ) - } - } - } - - return newCache[parameters][order]!! -} - /** * Compile the sizes array. * diff --git a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/expressions/DerivativeStructure.kt b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/expressions/DerivativeStructure.kt index a1a6354f0..01c045cdb 100644 --- a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/expressions/DerivativeStructure.kt +++ b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/expressions/DerivativeStructure.kt @@ -6,10 +6,9 @@ package space.kscience.kmath.expressions import space.kscience.kmath.misc.UnstableKMathAPI -import space.kscience.kmath.operations.NumericAlgebra import space.kscience.kmath.operations.Ring -import space.kscience.kmath.operations.ScaleOperations -import space.kscience.kmath.structures.MutableBuffer +import space.kscience.kmath.structures.Buffer +import space.kscience.kmath.structures.asBuffer /** * Class representing both the value and the differentials of a function. @@ -28,128 +27,29 @@ import space.kscience.kmath.structures.MutableBuffer * [Commons Math's `DerivativeStructure`](https://github.com/apache/commons-math/blob/924f6c357465b39beb50e3c916d5eb6662194175/commons-math-legacy/src/main/java/org/apache/commons/math4/legacy/analysis/differentiation/DerivativeStructure.java). */ @UnstableKMathAPI -public open class DerivativeStructure internal constructor( - internal val derivativeAlgebra: DerivativeStructureRing, - internal val compiler: DSCompiler, -) where A : Ring, A : NumericAlgebra, A : ScaleOperations { - /** - * Combined array holding all values. - */ - internal var data: MutableBuffer = - derivativeAlgebra.bufferFactory(compiler.size) { derivativeAlgebra.algebra.zero } +public open class DerivativeStructure> @PublishedApi internal constructor( + private val derivativeAlgebra: DerivativeStructureAlgebra, + @PublishedApi internal val data: Buffer, +) { - /** - * Build an instance with all values and derivatives set to 0. - * - * @param parameters number of free parameters. - * @param order derivation order. - */ - public constructor ( - derivativeAlgebra: DerivativeStructureRing, - parameters: Int, - order: Int, - ) : this( - derivativeAlgebra, - getCompiler(derivativeAlgebra.algebra, derivativeAlgebra.bufferFactory, parameters, order), - ) - - /** - * Build an instance representing a constant value. - * - * @param parameters number of free parameters. - * @param order derivation order. - * @param value value of the constant. - * @see DerivativeStructure - */ - public constructor ( - derivativeAlgebra: DerivativeStructureRing, - parameters: Int, - order: Int, - value: T, - ) : this( - derivativeAlgebra, - parameters, - order, - ) { - data[0] = value - } - - /** - * Build an instance representing a variable. - * - * Instances built using this constructor are considered to be the free variables with respect to which - * differentials are computed. As such, their differential with respect to themselves is +1. - * - * @param parameters number of free parameters. - * @param order derivation order. - * @param index index of the variable (from 0 to `parameters - 1`). - * @param value value of the variable. - */ - public constructor ( - derivativeAlgebra: DerivativeStructureRing, - parameters: Int, - order: Int, - index: Int, - value: T, - ) : this(derivativeAlgebra, parameters, order, value) { - require(index < parameters) { "number is too large: $index >= $parameters" } - - if (order > 0) { - // the derivative of the variable with respect to itself is 1. - data[getCompiler(derivativeAlgebra.algebra, derivativeAlgebra.bufferFactory, index, order).size] = - derivativeAlgebra.algebra.one - } - } - - /** - * Build an instance from all its derivatives. - * - * @param parameters number of free parameters. - * @param order derivation order. - * @param derivatives derivatives sorted according to [DSCompiler.getPartialDerivativeIndex]. - */ - public constructor ( - derivativeAlgebra: DerivativeStructureRing, - parameters: Int, - order: Int, - vararg derivatives: T, - ) : this( - derivativeAlgebra, - parameters, - order, - ) { - require(derivatives.size == data.size) { "dimension mismatch: ${derivatives.size} and ${data.size}" } - data = derivativeAlgebra.bufferFactory(data.size) { derivatives[it] } - } - - /** - * Copy constructor. - * - * @param ds instance to copy. - */ - internal constructor(ds: DerivativeStructure) : this(ds.derivativeAlgebra, ds.compiler) { - this.data = ds.data.copy() - } + public val compiler: DSCompiler get() = derivativeAlgebra.compiler /** * The number of free parameters. */ - public val freeParameters: Int - get() = compiler.freeParameters + public val freeParameters: Int get() = compiler.freeParameters /** * The derivation order. */ - public val order: Int - get() = compiler.order + public val order: Int get() = compiler.order /** * The value part of the derivative structure. * * @see getPartialDerivative */ - public val value: T - get() = data[0] + public val value: T get() = data[0] /** * Get a partial derivative. @@ -183,4 +83,75 @@ public open class DerivativeStructure internal constructor( public override fun hashCode(): Int = 227 + 229 * freeParameters + 233 * order + 239 * data.hashCode() + + public companion object { + + /** + * Build an instance representing a variable. + * + * Instances built using this constructor are considered to be the free variables with respect to which + * differentials are computed. As such, their differential with respect to themselves is +1. + */ + public fun > variable( + derivativeAlgebra: DerivativeStructureAlgebra, + index: Int, + value: T, + ): DerivativeStructure { + val compiler = derivativeAlgebra.compiler + require(index < compiler.freeParameters) { "number is too large: $index >= ${compiler.freeParameters}" } + return DerivativeStructure(derivativeAlgebra, derivativeAlgebra.bufferForVariable(index, value)) + } + + /** + * Build an instance from all its derivatives. + * + * @param derivatives derivatives sorted according to [DSCompiler.getPartialDerivativeIndex]. + */ + public fun > ofDerivatives( + derivativeAlgebra: DerivativeStructureAlgebra, + vararg derivatives: T, + ): DerivativeStructure { + val compiler = derivativeAlgebra.compiler + require(derivatives.size == compiler.size) { "dimension mismatch: ${derivatives.size} and ${compiler.size}" } + val data = derivatives.asBuffer() + + return DerivativeStructure( + derivativeAlgebra, + data + ) + } + } +} + +@OptIn(UnstableKMathAPI::class) +private fun > DerivativeStructureAlgebra.bufferForVariable(index: Int, value: T): Buffer { + val buffer = bufferFactory(compiler.size) { algebra.zero } + buffer[0] = value + if (compiler.order > 0) { + // the derivative of the variable with respect to itself is 1. + + val indexOfDerivative = compiler.getPartialDerivativeIndex(*IntArray(numberOfVariables).apply { + set(index, 1) + }) + + buffer[indexOfDerivative] = algebra.one + } + return buffer +} + +/** + * A class implementing both [DerivativeStructure] and [Symbol]. + */ +@UnstableKMathAPI +public class DerivativeStructureSymbol> internal constructor( + derivativeAlgebra: DerivativeStructureAlgebra, + index: Int, + symbol: Symbol, + value: T, +) : Symbol by symbol, DerivativeStructure( + derivativeAlgebra, derivativeAlgebra.bufferForVariable(index, value) +) { + override fun toString(): String = symbol.toString() + override fun equals(other: Any?): Boolean = (other as? Symbol) == symbol + override fun hashCode(): Int = symbol.hashCode() } diff --git a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/expressions/DerivativeStructureExpression.kt b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/expressions/DerivativeStructureExpression.kt index f91fb55e8..638057921 100644 --- a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/expressions/DerivativeStructureExpression.kt +++ b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/expressions/DerivativeStructureExpression.kt @@ -7,83 +7,89 @@ package space.kscience.kmath.expressions import space.kscience.kmath.misc.UnstableKMathAPI import space.kscience.kmath.operations.* +import space.kscience.kmath.structures.Buffer +import space.kscience.kmath.structures.MutableBuffer import space.kscience.kmath.structures.MutableBufferFactory -import space.kscience.kmath.structures.indices +import kotlin.math.max +import kotlin.math.min -/** - * A class implementing both [DerivativeStructure] and [Symbol]. - */ @UnstableKMathAPI -public class DerivativeStructureSymbol( - derivativeAlgebra: DerivativeStructureRing, - size: Int, - order: Int, - index: Int, - symbol: Symbol, - value: T, -) : Symbol by symbol, DerivativeStructure( - derivativeAlgebra, - size, - order, - index, - value -) where A : Ring, A : NumericAlgebra, A : ScaleOperations { - override fun toString(): String = symbol.toString() - override fun equals(other: Any?): Boolean = (other as? Symbol) == symbol - override fun hashCode(): Int = symbol.hashCode() -} - -/** - * A ring over [DerivativeStructure]. - * - * @property order The derivation order. - * @param bindings The map of bindings values. All bindings are considered free parameters. - */ -@UnstableKMathAPI -public open class DerivativeStructureRing( +public abstract class DerivativeStructureAlgebra>( public val algebra: A, public val bufferFactory: MutableBufferFactory, public val order: Int, bindings: Map, -) : Ring>, ScaleOperations>, - NumericAlgebra>, - ExpressionAlgebra>, - NumbersAddOps> where A : Ring, A : NumericAlgebra, A : ScaleOperations { +) : ExpressionAlgebra> { + public val numberOfVariables: Int = bindings.size - override val zero: DerivativeStructure by lazy { - DerivativeStructure( - this, - numberOfVariables, - order, - ) - } - override val one: DerivativeStructure by lazy { - DerivativeStructure( - this, - numberOfVariables, - order, - algebra.one, - ) - } + /** + * Get the compiler for number of free parameters and order. + * + * @return cached rules set. + */ + @PublishedApi + internal val compiler: DSCompiler by lazy { + // get the cached compilers + val cache: Array?>>? = null - override fun number(value: Number): DerivativeStructure = const(algebra.number(value)) + // we need to create more compilers + val maxParameters: Int = max(numberOfVariables, cache?.size ?: 0) + val maxOrder: Int = max(order, if (cache == null) 0 else cache[0].size) + val newCache: Array?>> = Array(maxParameters + 1) { arrayOfNulls(maxOrder + 1) } + + if (cache != null) { + // preserve the already created compilers + for (i in cache.indices) { + cache[i].copyInto(newCache[i], endIndex = cache[i].size) + } + } + + // create the array in increasing diagonal order + for (diag in 0..numberOfVariables + order) { + for (o in max(0, diag - numberOfVariables)..min(order, diag)) { + val p: Int = diag - o + if (newCache[p][o] == null) { + val valueCompiler: DSCompiler? = if (p == 0) null else newCache[p - 1][o]!! + val derivativeCompiler: DSCompiler? = if (o == 0) null else newCache[p][o - 1]!! + + newCache[p][o] = DSCompiler( + algebra, + bufferFactory, + p, + o, + valueCompiler, + derivativeCompiler, + ) + } + } + } + + return@lazy newCache[numberOfVariables][order]!! + } private val variables: Map> = bindings.entries.mapIndexed { index, (key, value) -> key to DerivativeStructureSymbol( this, - numberOfVariables, - order, index, key, value, ) }.toMap() - public override fun const(value: T): DerivativeStructure = - DerivativeStructure(this, numberOfVariables, order, value) + + + public override fun const(value: T): DerivativeStructure { + val buffer = bufferFactory(compiler.size) { algebra.zero } + buffer[0] = value + + return DerivativeStructure( + this, + buffer + ) + } override fun bindSymbolOrNull(value: String): DerivativeStructureSymbol? = variables[StringSymbol(value)] @@ -103,54 +109,99 @@ public open class DerivativeStructureRing( public fun DerivativeStructure.derivative(vararg symbols: Symbol): T = derivative(symbols.toList()) +} + + +/** + * A ring over [DerivativeStructure]. + * + * @property order The derivation order. + * @param bindings The map of bindings values. All bindings are considered free parameters. + */ +@UnstableKMathAPI +public open class DerivativeStructureRing( + algebra: A, + bufferFactory: MutableBufferFactory, + order: Int, + bindings: Map, +) : DerivativeStructureAlgebra(algebra, bufferFactory, order, bindings), + Ring>, ScaleOperations>, + NumericAlgebra>, + NumbersAddOps> where A : Ring, A : NumericAlgebra, A : ScaleOperations { + + override fun bindSymbolOrNull(value: String): DerivativeStructureSymbol? = + super.bindSymbolOrNull(value) + override fun DerivativeStructure.unaryMinus(): DerivativeStructure { - val ds = DerivativeStructure(this@DerivativeStructureRing, compiler) - for (i in ds.data.indices) { - ds.data[i] = algebra { -data[i] } - } - return ds + val newData = algebra { data.map(bufferFactory) { -it } } + return DerivativeStructure(this@DerivativeStructureRing, newData) } + /** + * Create a copy of given [Buffer] and modify it according to [block] + */ + protected inline fun DerivativeStructure.transformDataBuffer(block: DSCompiler.(MutableBuffer) -> Unit): DerivativeStructure { + val newData = bufferFactory(compiler.size) { data[it] } + compiler.block(newData) + return DerivativeStructure(this@DerivativeStructureRing, newData) + } + + protected fun DerivativeStructure.mapData(block: (T) -> T): DerivativeStructure { + val newData: Buffer = data.map(bufferFactory, block) + return DerivativeStructure(this@DerivativeStructureRing, newData) + } + + protected fun DerivativeStructure.mapDataIndexed(block: (Int, T) -> T): DerivativeStructure { + val newData: Buffer = data.mapIndexed(bufferFactory, block) + return DerivativeStructure(this@DerivativeStructureRing, newData) + } + + override val zero: DerivativeStructure by lazy { + const(algebra.zero) + } + + override val one: DerivativeStructure by lazy { + const(algebra.one) + } + + override fun number(value: Number): DerivativeStructure = const(algebra.number(value)) + override fun add(left: DerivativeStructure, right: DerivativeStructure): DerivativeStructure { left.compiler.checkCompatibility(right.compiler) - val ds = DerivativeStructure(left) - left.compiler.add(left.data, 0, right.data, 0, ds.data, 0) - return ds + return left.transformDataBuffer { result -> + add(left.data, 0, right.data, 0, result, 0) + } } - override fun scale(a: DerivativeStructure, value: Double): DerivativeStructure { - val ds = DerivativeStructure(a) - for (i in ds.data.indices) { - ds.data[i] = algebra { ds.data[i].times(value) } - } - return ds + override fun scale(a: DerivativeStructure, value: Double): DerivativeStructure = algebra { + a.mapData { it.times(value) } } override fun multiply( left: DerivativeStructure, - right: DerivativeStructure + right: DerivativeStructure, ): DerivativeStructure { left.compiler.checkCompatibility(right.compiler) - val result = DerivativeStructure(this, left.compiler) - left.compiler.multiply(left.data, 0, right.data, 0, result.data, 0) - return result + return left.transformDataBuffer { result -> + multiply(left.data, 0, right.data, 0, result, 0) + } } override fun DerivativeStructure.minus(arg: DerivativeStructure): DerivativeStructure { compiler.checkCompatibility(arg.compiler) - val ds = DerivativeStructure(this) - compiler.subtract(data, 0, arg.data, 0, ds.data, 0) - return ds + return transformDataBuffer { result -> + subtract(data, 0, arg.data, 0, result, 0) + } } - override operator fun DerivativeStructure.plus(other: Number): DerivativeStructure { - val ds = DerivativeStructure(this) - ds.data[0] = algebra { ds.data[0] + number(other) } - return ds + override operator fun DerivativeStructure.plus(other: Number): DerivativeStructure = algebra { + transformDataBuffer { + it[0] += number(other) + } } override operator fun DerivativeStructure.minus(other: Number): DerivativeStructure = - this + -other.toDouble() + this + (-other.toDouble()) override operator fun Number.plus(other: DerivativeStructure): DerivativeStructure = other + this override operator fun Number.minus(other: DerivativeStructure): DerivativeStructure = other - this @@ -194,119 +245,85 @@ public class DerivativeStructureField>( override fun divide(left: DerivativeStructure, right: DerivativeStructure): DerivativeStructure { left.compiler.checkCompatibility(right.compiler) - val result = DerivativeStructure(this, left.compiler) - left.compiler.divide(left.data, 0, right.data, 0, result.data, 0) - return result + return left.transformDataBuffer { result -> + left.compiler.divide(left.data, 0, right.data, 0, result, 0) + } } - override fun sin(arg: DerivativeStructure): DerivativeStructure { - val result = DerivativeStructure(this, arg.compiler) - arg.compiler.sin(arg.data, 0, result.data, 0) - return result + override fun sin(arg: DerivativeStructure): DerivativeStructure = arg.transformDataBuffer { result -> + sin(arg.data, 0, result, 0) } - override fun cos(arg: DerivativeStructure): DerivativeStructure { - val result = DerivativeStructure(this, arg.compiler) - arg.compiler.cos(arg.data, 0, result.data, 0) - return result + override fun cos(arg: DerivativeStructure): DerivativeStructure = arg.transformDataBuffer { result -> + cos(arg.data, 0, result, 0) } - override fun tan(arg: DerivativeStructure): DerivativeStructure { - val result = DerivativeStructure(this, arg.compiler) - arg.compiler.tan(arg.data, 0, result.data, 0) - return result + override fun tan(arg: DerivativeStructure): DerivativeStructure = arg.transformDataBuffer { result -> + tan(arg.data, 0, result, 0) } - override fun asin(arg: DerivativeStructure): DerivativeStructure { - val result = DerivativeStructure(this, arg.compiler) - arg.compiler.asin(arg.data, 0, result.data, 0) - return result + override fun asin(arg: DerivativeStructure): DerivativeStructure = arg.transformDataBuffer { result -> + asin(arg.data, 0, result, 0) } - override fun acos(arg: DerivativeStructure): DerivativeStructure { - val result = DerivativeStructure(this, arg.compiler) - arg.compiler.acos(arg.data, 0, result.data, 0) - return result + override fun acos(arg: DerivativeStructure): DerivativeStructure = arg.transformDataBuffer { result -> + acos(arg.data, 0, result, 0) } - override fun atan(arg: DerivativeStructure): DerivativeStructure { - val result = DerivativeStructure(this, arg.compiler) - arg.compiler.atan(arg.data, 0, result.data, 0) - return result + override fun atan(arg: DerivativeStructure): DerivativeStructure = arg.transformDataBuffer { result -> + atan(arg.data, 0, result, 0) } - override fun sinh(arg: DerivativeStructure): DerivativeStructure { - val result = DerivativeStructure(this, arg.compiler) - arg.compiler.sinh(arg.data, 0, result.data, 0) - return result + override fun sinh(arg: DerivativeStructure): DerivativeStructure = arg.transformDataBuffer { result -> + sinh(arg.data, 0, result, 0) } - override fun cosh(arg: DerivativeStructure): DerivativeStructure { - val result = DerivativeStructure(this, arg.compiler) - arg.compiler.cosh(arg.data, 0, result.data, 0) - return result + override fun cosh(arg: DerivativeStructure): DerivativeStructure = arg.transformDataBuffer { result -> + cosh(arg.data, 0, result, 0) } - override fun tanh(arg: DerivativeStructure): DerivativeStructure { - val result = DerivativeStructure(this, arg.compiler) - arg.compiler.tanh(arg.data, 0, result.data, 0) - return result + override fun tanh(arg: DerivativeStructure): DerivativeStructure = arg.transformDataBuffer { result -> + tanh(arg.data, 0, result, 0) } - override fun asinh(arg: DerivativeStructure): DerivativeStructure { - val result = DerivativeStructure(this, arg.compiler) - arg.compiler.asinh(arg.data, 0, result.data, 0) - return result + override fun asinh(arg: DerivativeStructure): DerivativeStructure = arg.transformDataBuffer { result -> + asinh(arg.data, 0, result, 0) } - override fun acosh(arg: DerivativeStructure): DerivativeStructure { - val result = DerivativeStructure(this, arg.compiler) - arg.compiler.acosh(arg.data, 0, result.data, 0) - return result + override fun acosh(arg: DerivativeStructure): DerivativeStructure = arg.transformDataBuffer { result -> + acosh(arg.data, 0, result, 0) } - override fun atanh(arg: DerivativeStructure): DerivativeStructure { - val result = DerivativeStructure(this, arg.compiler) - arg.compiler.atanh(arg.data, 0, result.data, 0) - return result + override fun atanh(arg: DerivativeStructure): DerivativeStructure = arg.transformDataBuffer { result -> + atanh(arg.data, 0, result, 0) } override fun power(arg: DerivativeStructure, pow: Number): DerivativeStructure = when (pow) { - is Int -> { - val result = DerivativeStructure(this, arg.compiler) - arg.compiler.pow(arg.data, 0, pow, result.data, 0) - result + is Int -> arg.transformDataBuffer { result -> + pow(arg.data, 0, pow, result, 0) } - else -> { - val result = DerivativeStructure(this, arg.compiler) - arg.compiler.pow(arg.data, 0, pow.toDouble(), result.data, 0) - result + else -> arg.transformDataBuffer { result -> + pow(arg.data, 0, pow.toDouble(), result, 0) } } - override fun sqrt(arg: DerivativeStructure): DerivativeStructure { - val result = DerivativeStructure(this, arg.compiler) - arg.compiler.sqrt(arg.data, 0, result.data, 0) - return result + override fun sqrt(arg: DerivativeStructure): DerivativeStructure = arg.transformDataBuffer { result -> + sqrt(arg.data, 0, result, 0) } public fun power(arg: DerivativeStructure, pow: DerivativeStructure): DerivativeStructure { arg.compiler.checkCompatibility(pow.compiler) - val result = DerivativeStructure(this, arg.compiler) - arg.compiler.pow(arg.data, 0, pow.data, 0, result.data, 0) - return result + return arg.transformDataBuffer { result -> + pow(arg.data, 0, pow.data, 0, result, 0) + } } - override fun exp(arg: DerivativeStructure): DerivativeStructure { - val result = DerivativeStructure(this, arg.compiler) - arg.compiler.exp(arg.data, 0, result.data, 0) - return result + override fun exp(arg: DerivativeStructure): DerivativeStructure = arg.transformDataBuffer { result -> + exp(arg.data, 0, result, 0) } - override fun ln(arg: DerivativeStructure): DerivativeStructure { - val result = DerivativeStructure(this, arg.compiler) - arg.compiler.ln(arg.data, 0, result.data, 0) - return result + override fun ln(arg: DerivativeStructure): DerivativeStructure = arg.transformDataBuffer { result -> + ln(arg.data, 0, result, 0) } } 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 715fad07b..10438dd02 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 @@ -188,7 +188,7 @@ public interface LinearSpace> { */ public fun > buffered( algebra: A, - bufferFactory: BufferFactory = Buffer.Companion::boxing, + bufferFactory: BufferFactory = BufferFactory(Buffer.Companion::boxing), ): LinearSpace = BufferedLinearSpace(BufferRingOps(algebra, bufferFactory)) @Deprecated("use DoubleField.linearSpace") diff --git a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/misc/annotations.kt b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/misc/annotations.kt index 7c612b6a9..60fa81cd8 100644 --- a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/misc/annotations.kt +++ b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/misc/annotations.kt @@ -27,5 +27,5 @@ public annotation class UnstableKMathAPI RequiresOptIn.Level.WARNING, ) public annotation class PerformancePitfall( - val message: String = "Potential performance problem" + val message: String = "Potential performance problem", ) diff --git a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/BufferND.kt b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/BufferND.kt index 2401f6319..8175bd65e 100644 --- a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/BufferND.kt +++ b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/BufferND.kt @@ -69,7 +69,7 @@ public class MutableBufferND( * Transform structure to a new structure using provided [MutableBufferFactory] and optimizing if argument is [MutableBufferND] */ public inline fun MutableStructureND.mapToMutableBuffer( - factory: MutableBufferFactory = MutableBuffer.Companion::auto, + factory: MutableBufferFactory = MutableBufferFactory(MutableBuffer.Companion::auto), crossinline transform: (T) -> R, ): MutableBufferND { return if (this is MutableBufferND) diff --git a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/StructureND.kt b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/StructureND.kt index e934c6370..6e54e1b9d 100644 --- a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/StructureND.kt +++ b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/StructureND.kt @@ -120,7 +120,7 @@ public interface StructureND : Featured, WithShape { */ public fun buffered( strides: Strides, - bufferFactory: BufferFactory = Buffer.Companion::boxing, + bufferFactory: BufferFactory = BufferFactory(Buffer.Companion::boxing), initializer: (IntArray) -> T, ): BufferND = BufferND(strides, bufferFactory(strides.linearSize) { i -> initializer(strides.index(i)) }) @@ -140,7 +140,7 @@ public interface StructureND : Featured, WithShape { public fun buffered( shape: IntArray, - bufferFactory: BufferFactory = Buffer.Companion::boxing, + bufferFactory: BufferFactory = BufferFactory(Buffer.Companion::boxing), initializer: (IntArray) -> T, ): BufferND = buffered(DefaultStrides(shape), bufferFactory, initializer) diff --git a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/operations/DoubleBufferOps.kt b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/operations/DoubleBufferOps.kt index 0ee591acc..083892105 100644 --- a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/operations/DoubleBufferOps.kt +++ b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/operations/DoubleBufferOps.kt @@ -6,12 +6,10 @@ package space.kscience.kmath.operations import space.kscience.kmath.linear.Point -import space.kscience.kmath.misc.UnstableKMathAPI import space.kscience.kmath.structures.Buffer import space.kscience.kmath.structures.BufferFactory import space.kscience.kmath.structures.DoubleBuffer import space.kscience.kmath.structures.asBuffer - import kotlin.math.* /** @@ -21,7 +19,7 @@ public abstract class DoubleBufferOps : BufferAlgebra, Exte Norm, Double> { override val elementAlgebra: DoubleField get() = DoubleField - override val bufferFactory: BufferFactory get() = ::DoubleBuffer + override val bufferFactory: BufferFactory get() = BufferFactory(::DoubleBuffer) override fun Buffer.map(block: DoubleField.(Double) -> Double): DoubleBuffer = mapInline { DoubleField.block(it) } diff --git a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/operations/bufferOperation.kt b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/operations/bufferOperation.kt index 31b0c2841..652472fcf 100644 --- a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/operations/bufferOperation.kt +++ b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/operations/bufferOperation.kt @@ -61,31 +61,39 @@ public inline fun Buffer.toTypedArray(): Array = Array(size, : /** * Create a new buffer from this one with the given mapping function and using [Buffer.Companion.auto] buffer factory. */ -public inline fun Buffer.map(block: (T) -> R): Buffer = +public inline fun Buffer.map(block: (T) -> R): Buffer = Buffer.auto(size) { block(get(it)) } /** * Create a new buffer from this one with the given mapping function. * Provided [bufferFactory] is used to construct the new buffer. */ -public inline fun Buffer.map( +public inline fun Buffer.map( bufferFactory: BufferFactory, crossinline block: (T) -> R, ): Buffer = bufferFactory(size) { block(get(it)) } /** - * Create a new buffer from this one with the given indexed mapping function. - * Provided [BufferFactory] is used to construct the new buffer. + * Create a new buffer from this one with the given mapping (indexed) function. + * Provided [bufferFactory] is used to construct the new buffer. */ -public inline fun Buffer.mapIndexed( - bufferFactory: BufferFactory = Buffer.Companion::auto, +public inline fun Buffer.mapIndexed( + bufferFactory: BufferFactory, crossinline block: (index: Int, value: T) -> R, ): Buffer = bufferFactory(size) { block(it, get(it)) } +/** + * Create a new buffer from this one with the given indexed mapping function. + * Provided [BufferFactory] is used to construct the new buffer. + */ +public inline fun Buffer.mapIndexed( + crossinline block: (index: Int, value: T) -> R, +): Buffer = BufferFactory(Buffer.Companion::auto).invoke(size) { block(it, get(it)) } + /** * Fold given buffer according to [operation] */ -public inline fun Buffer.fold(initial: R, operation: (acc: R, T) -> R): R { +public inline fun Buffer.fold(initial: R, operation: (acc: R, T) -> R): R { var accumulator = initial for (index in this.indices) accumulator = operation(accumulator, get(index)) return accumulator @@ -95,9 +103,9 @@ public inline fun Buffer.fold(initial: R, operation: (acc: R, T) * Zip two buffers using given [transform]. */ @UnstableKMathAPI -public inline fun Buffer.zip( +public inline fun Buffer.zip( other: Buffer, - bufferFactory: BufferFactory = Buffer.Companion::auto, + bufferFactory: BufferFactory = BufferFactory(Buffer.Companion::auto), crossinline transform: (T1, T2) -> R, ): Buffer { require(size == other.size) { "Buffer size mismatch in zip: expected $size but found ${other.size}" } diff --git a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/structures/Buffer.kt b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/structures/Buffer.kt index a1b0307c4..1c79c257a 100644 --- a/kmath-core/src/commonMain/kotlin/space/kscience/kmath/structures/Buffer.kt +++ b/kmath-core/src/commonMain/kotlin/space/kscience/kmath/structures/Buffer.kt @@ -14,14 +14,18 @@ import kotlin.reflect.KClass * * @param T the type of buffer. */ -public typealias BufferFactory = (Int, (Int) -> T) -> Buffer +public fun interface BufferFactory { + public operator fun invoke(size: Int, builder: (Int) -> T): Buffer +} /** * Function that produces [MutableBuffer] from its size and function that supplies values. * * @param T the type of buffer. */ -public typealias MutableBufferFactory = (Int, (Int) -> T) -> MutableBuffer +public fun interface MutableBufferFactory: BufferFactory{ + override fun invoke(size: Int, builder: (Int) -> T): MutableBuffer +} /** * A generic read-only random-access structure for both primitives and objects. diff --git a/kmath-ejml/src/main/kotlin/space/kscience/kmath/ejml/_generated.kt b/kmath-ejml/src/main/kotlin/space/kscience/kmath/ejml/_generated.kt deleted file mode 100644 index aac327a84..000000000 --- a/kmath-ejml/src/main/kotlin/space/kscience/kmath/ejml/_generated.kt +++ /dev/null @@ -1,1003 +0,0 @@ -/* - * 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 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(override val origin: M) : EjmlVector(origin) { - init { - require(origin.numRows == 1) { "The origin matrix must have only one row to form a vector" } - } - - override operator fun get(index: Int): Double = origin[0, index] -} - -/** - * [EjmlVector] specialization for [Float]. - */ -public class EjmlFloatVector(override val origin: M) : EjmlVector(origin) { - init { - require(origin.numRows == 1) { "The origin matrix must have only one row to form a vector" } - } - - override operator fun get(index: Int): Float = origin[0, index] -} - -/** - * [EjmlMatrix] specialization for [Double]. - */ -public class EjmlDoubleMatrix(override val origin: M) : EjmlMatrix(origin) { - override operator fun get(i: Int, j: Int): Double = origin[i, j] -} - -/** - * [EjmlMatrix] specialization for [Float]. - */ -public class EjmlFloatMatrix(override val origin: M) : EjmlMatrix(origin) { - 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() { - /** - * The [DoubleField] reference. - */ - override val elementAlgebra: DoubleField get() = DoubleField - - @Suppress("UNCHECKED_CAST") - override fun Matrix.toEjml(): EjmlDoubleMatrix = when { - this is EjmlDoubleMatrix<*> && origin is DMatrixRMaj -> this as EjmlDoubleMatrix - else -> buildMatrix(rowNum, colNum) { i, j -> get(i, j) } - } - - @Suppress("UNCHECKED_CAST") - override fun Point.toEjml(): EjmlDoubleVector = when { - this is EjmlDoubleVector<*> && origin is DMatrixRMaj -> this as EjmlDoubleVector - else -> EjmlDoubleVector(DMatrixRMaj(size, 1).also { - (0 until it.numRows).forEach { row -> it[row, 0] = get(row) } - }) - } - - override fun buildMatrix( - rows: Int, - columns: Int, - initializer: DoubleField.(i: Int, j: Int) -> Double, - ): EjmlDoubleMatrix = DMatrixRMaj(rows, columns).also { - (0 until rows).forEach { row -> - (0 until columns).forEach { col -> it[row, col] = elementAlgebra.initializer(row, col) } - } - }.wrapMatrix() - - override fun buildVector( - size: Int, - initializer: DoubleField.(Int) -> Double, - ): EjmlDoubleVector = EjmlDoubleVector(DMatrixRMaj(size, 1).also { - (0 until it.numRows).forEach { row -> it[row, 0] = elementAlgebra.initializer(row) } - }) - - private fun T.wrapMatrix() = EjmlDoubleMatrix(this) - private fun T.wrapVector() = EjmlDoubleVector(this) - - override fun Matrix.unaryMinus(): Matrix = this * elementAlgebra { -one } - - override fun Matrix.dot(other: Matrix): EjmlDoubleMatrix { - val out = DMatrixRMaj(1, 1) - CommonOps_DDRM.mult(toEjml().origin, other.toEjml().origin, out) - return out.wrapMatrix() - } - - override fun Matrix.dot(vector: Point): EjmlDoubleVector { - val out = DMatrixRMaj(1, 1) - CommonOps_DDRM.mult(toEjml().origin, vector.toEjml().origin, out) - return out.wrapVector() - } - - override operator fun Matrix.minus(other: Matrix): EjmlDoubleMatrix { - val out = DMatrixRMaj(1, 1) - - CommonOps_DDRM.add( - elementAlgebra.one, - toEjml().origin, - elementAlgebra { -one }, - other.toEjml().origin, - out, - ) - - return out.wrapMatrix() - } - - override operator fun Matrix.times(value: Double): EjmlDoubleMatrix { - val res = DMatrixRMaj(1, 1) - CommonOps_DDRM.scale(value, toEjml().origin, res) - return res.wrapMatrix() - } - - override fun Point.unaryMinus(): EjmlDoubleVector { - val res = DMatrixRMaj(1, 1) - CommonOps_DDRM.changeSign(toEjml().origin, res) - return res.wrapVector() - } - - override fun Matrix.plus(other: Matrix): EjmlDoubleMatrix { - val out = DMatrixRMaj(1, 1) - - CommonOps_DDRM.add( - elementAlgebra.one, - toEjml().origin, - elementAlgebra.one, - other.toEjml().origin, - out, - ) - - return out.wrapMatrix() - } - - override fun Point.plus(other: Point): EjmlDoubleVector { - val out = DMatrixRMaj(1, 1) - - CommonOps_DDRM.add( - elementAlgebra.one, - toEjml().origin, - elementAlgebra.one, - other.toEjml().origin, - out, - ) - - return out.wrapVector() - } - - override fun Point.minus(other: Point): EjmlDoubleVector { - val out = DMatrixRMaj(1, 1) - - CommonOps_DDRM.add( - elementAlgebra.one, - toEjml().origin, - elementAlgebra { -one }, - other.toEjml().origin, - out, - ) - - return out.wrapVector() - } - - override fun Double.times(m: Matrix): EjmlDoubleMatrix = m * this - - override fun Point.times(value: Double): EjmlDoubleVector { - val res = DMatrixRMaj(1, 1) - CommonOps_DDRM.scale(value, toEjml().origin, res) - return res.wrapVector() - } - - override fun Double.times(v: Point): EjmlDoubleVector = v * this - - @UnstableKMathAPI - override fun computeFeature(structure: Matrix, type: KClass): F? { - structure.getFeature(type)?.let { return it } - val origin = structure.toEjml().origin - - return when (type) { - InverseMatrixFeature::class -> object : InverseMatrixFeature { - override val inverse: Matrix by lazy { - val res = origin.copy() - CommonOps_DDRM.invert(res) - res.wrapMatrix() - } - } - - DeterminantFeature::class -> object : DeterminantFeature { - override val determinant: Double by lazy { CommonOps_DDRM.det(origin) } - } - - SingularValueDecompositionFeature::class -> object : SingularValueDecompositionFeature { - private val svd by lazy { - DecompositionFactory_DDRM.svd(origin.numRows, origin.numCols, true, true, false) - .apply { decompose(origin.copy()) } - } - - override val u: Matrix by lazy { svd.getU(null, false).wrapMatrix() } - override val s: Matrix by lazy { svd.getW(null).wrapMatrix() } - override val v: Matrix by lazy { svd.getV(null, false).wrapMatrix() } - override val singularValues: Point by lazy { DoubleBuffer(svd.singularValues) } - } - - QRDecompositionFeature::class -> object : QRDecompositionFeature { - private val qr by lazy { - DecompositionFactory_DDRM.qr().apply { decompose(origin.copy()) } - } - - override val q: Matrix by lazy { - qr.getQ(null, false).wrapMatrix().withFeature(OrthogonalFeature) - } - - override val r: Matrix by lazy { qr.getR(null, false).wrapMatrix().withFeature(UFeature) } - } - - CholeskyDecompositionFeature::class -> object : CholeskyDecompositionFeature { - override val l: Matrix by lazy { - val cholesky = - DecompositionFactory_DDRM.chol(structure.rowNum, true).apply { decompose(origin.copy()) } - - cholesky.getT(null).wrapMatrix().withFeature(LFeature) - } - } - - LupDecompositionFeature::class -> object : LupDecompositionFeature { - private val lup by lazy { - DecompositionFactory_DDRM.lu(origin.numRows, origin.numCols).apply { decompose(origin.copy()) } - } - - override val l: Matrix by lazy { - lup.getLower(null).wrapMatrix().withFeature(LFeature) - } - - override val u: Matrix by lazy { - lup.getUpper(null).wrapMatrix().withFeature(UFeature) - } - - override val p: Matrix by lazy { lup.getRowPivot(null).wrapMatrix() } - } - - else -> null - }?.let{ - type.cast(it) - } - } - - /** - * Solves for *x* in the following equation: *x = [a] -1 · [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, b: Matrix): EjmlDoubleMatrix { - 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] -1 · [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, b: Point): EjmlDoubleVector { - 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() { - /** - * The [FloatField] reference. - */ - override val elementAlgebra: FloatField get() = FloatField - - @Suppress("UNCHECKED_CAST") - override fun Matrix.toEjml(): EjmlFloatMatrix = when { - this is EjmlFloatMatrix<*> && origin is FMatrixRMaj -> this as EjmlFloatMatrix - else -> buildMatrix(rowNum, colNum) { i, j -> get(i, j) } - } - - @Suppress("UNCHECKED_CAST") - override fun Point.toEjml(): EjmlFloatVector = when { - this is EjmlFloatVector<*> && origin is FMatrixRMaj -> this as EjmlFloatVector - else -> EjmlFloatVector(FMatrixRMaj(size, 1).also { - (0 until it.numRows).forEach { row -> it[row, 0] = get(row) } - }) - } - - override fun buildMatrix( - rows: Int, - columns: Int, - initializer: FloatField.(i: Int, j: Int) -> Float, - ): EjmlFloatMatrix = FMatrixRMaj(rows, columns).also { - (0 until rows).forEach { row -> - (0 until columns).forEach { col -> it[row, col] = elementAlgebra.initializer(row, col) } - } - }.wrapMatrix() - - override fun buildVector( - size: Int, - initializer: FloatField.(Int) -> Float, - ): EjmlFloatVector = EjmlFloatVector(FMatrixRMaj(size, 1).also { - (0 until it.numRows).forEach { row -> it[row, 0] = elementAlgebra.initializer(row) } - }) - - private fun T.wrapMatrix() = EjmlFloatMatrix(this) - private fun T.wrapVector() = EjmlFloatVector(this) - - override fun Matrix.unaryMinus(): Matrix = this * elementAlgebra { -one } - - override fun Matrix.dot(other: Matrix): EjmlFloatMatrix { - val out = FMatrixRMaj(1, 1) - CommonOps_FDRM.mult(toEjml().origin, other.toEjml().origin, out) - return out.wrapMatrix() - } - - override fun Matrix.dot(vector: Point): EjmlFloatVector { - val out = FMatrixRMaj(1, 1) - CommonOps_FDRM.mult(toEjml().origin, vector.toEjml().origin, out) - return out.wrapVector() - } - - override operator fun Matrix.minus(other: Matrix): EjmlFloatMatrix { - val out = FMatrixRMaj(1, 1) - - CommonOps_FDRM.add( - elementAlgebra.one, - toEjml().origin, - elementAlgebra { -one }, - other.toEjml().origin, - out, - ) - - return out.wrapMatrix() - } - - override operator fun Matrix.times(value: Float): EjmlFloatMatrix { - val res = FMatrixRMaj(1, 1) - CommonOps_FDRM.scale(value, toEjml().origin, res) - return res.wrapMatrix() - } - - override fun Point.unaryMinus(): EjmlFloatVector { - val res = FMatrixRMaj(1, 1) - CommonOps_FDRM.changeSign(toEjml().origin, res) - return res.wrapVector() - } - - override fun Matrix.plus(other: Matrix): EjmlFloatMatrix { - val out = FMatrixRMaj(1, 1) - - CommonOps_FDRM.add( - elementAlgebra.one, - toEjml().origin, - elementAlgebra.one, - other.toEjml().origin, - out, - ) - - return out.wrapMatrix() - } - - override fun Point.plus(other: Point): EjmlFloatVector { - val out = FMatrixRMaj(1, 1) - - CommonOps_FDRM.add( - elementAlgebra.one, - toEjml().origin, - elementAlgebra.one, - other.toEjml().origin, - out, - ) - - return out.wrapVector() - } - - override fun Point.minus(other: Point): EjmlFloatVector { - val out = FMatrixRMaj(1, 1) - - CommonOps_FDRM.add( - elementAlgebra.one, - toEjml().origin, - elementAlgebra { -one }, - other.toEjml().origin, - out, - ) - - return out.wrapVector() - } - - override fun Float.times(m: Matrix): EjmlFloatMatrix = m * this - - override fun Point.times(value: Float): EjmlFloatVector { - val res = FMatrixRMaj(1, 1) - CommonOps_FDRM.scale(value, toEjml().origin, res) - return res.wrapVector() - } - - override fun Float.times(v: Point): EjmlFloatVector = v * this - - @UnstableKMathAPI - override fun computeFeature(structure: Matrix, type: KClass): F? { - structure.getFeature(type)?.let { return it } - val origin = structure.toEjml().origin - - return when (type) { - InverseMatrixFeature::class -> object : InverseMatrixFeature { - override val inverse: Matrix by lazy { - val res = origin.copy() - CommonOps_FDRM.invert(res) - res.wrapMatrix() - } - } - - DeterminantFeature::class -> object : DeterminantFeature { - override val determinant: Float by lazy { CommonOps_FDRM.det(origin) } - } - - SingularValueDecompositionFeature::class -> object : SingularValueDecompositionFeature { - private val svd by lazy { - DecompositionFactory_FDRM.svd(origin.numRows, origin.numCols, true, true, false) - .apply { decompose(origin.copy()) } - } - - override val u: Matrix by lazy { svd.getU(null, false).wrapMatrix() } - override val s: Matrix by lazy { svd.getW(null).wrapMatrix() } - override val v: Matrix by lazy { svd.getV(null, false).wrapMatrix() } - override val singularValues: Point by lazy { FloatBuffer(svd.singularValues) } - } - - QRDecompositionFeature::class -> object : QRDecompositionFeature { - private val qr by lazy { - DecompositionFactory_FDRM.qr().apply { decompose(origin.copy()) } - } - - override val q: Matrix by lazy { - qr.getQ(null, false).wrapMatrix().withFeature(OrthogonalFeature) - } - - override val r: Matrix by lazy { qr.getR(null, false).wrapMatrix().withFeature(UFeature) } - } - - CholeskyDecompositionFeature::class -> object : CholeskyDecompositionFeature { - override val l: Matrix by lazy { - val cholesky = - DecompositionFactory_FDRM.chol(structure.rowNum, true).apply { decompose(origin.copy()) } - - cholesky.getT(null).wrapMatrix().withFeature(LFeature) - } - } - - LupDecompositionFeature::class -> object : LupDecompositionFeature { - private val lup by lazy { - DecompositionFactory_FDRM.lu(origin.numRows, origin.numCols).apply { decompose(origin.copy()) } - } - - override val l: Matrix by lazy { - lup.getLower(null).wrapMatrix().withFeature(LFeature) - } - - override val u: Matrix by lazy { - lup.getUpper(null).wrapMatrix().withFeature(UFeature) - } - - override val p: Matrix by lazy { lup.getRowPivot(null).wrapMatrix() } - } - - else -> null - }?.let{ - type.cast(it) - } - } - - /** - * Solves for *x* in the following equation: *x = [a] -1 · [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, b: Matrix): EjmlFloatMatrix { - 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] -1 · [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, b: Point): EjmlFloatVector { - 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() { - /** - * The [DoubleField] reference. - */ - override val elementAlgebra: DoubleField get() = DoubleField - - @Suppress("UNCHECKED_CAST") - override fun Matrix.toEjml(): EjmlDoubleMatrix = when { - this is EjmlDoubleMatrix<*> && origin is DMatrixSparseCSC -> this as EjmlDoubleMatrix - else -> buildMatrix(rowNum, colNum) { i, j -> get(i, j) } - } - - @Suppress("UNCHECKED_CAST") - override fun Point.toEjml(): EjmlDoubleVector = when { - this is EjmlDoubleVector<*> && origin is DMatrixSparseCSC -> this as EjmlDoubleVector - else -> EjmlDoubleVector(DMatrixSparseCSC(size, 1).also { - (0 until it.numRows).forEach { row -> it[row, 0] = get(row) } - }) - } - - override fun buildMatrix( - rows: Int, - columns: Int, - initializer: DoubleField.(i: Int, j: Int) -> Double, - ): EjmlDoubleMatrix = DMatrixSparseCSC(rows, columns).also { - (0 until rows).forEach { row -> - (0 until columns).forEach { col -> it[row, col] = elementAlgebra.initializer(row, col) } - } - }.wrapMatrix() - - override fun buildVector( - size: Int, - initializer: DoubleField.(Int) -> Double, - ): EjmlDoubleVector = EjmlDoubleVector(DMatrixSparseCSC(size, 1).also { - (0 until it.numRows).forEach { row -> it[row, 0] = elementAlgebra.initializer(row) } - }) - - private fun T.wrapMatrix() = EjmlDoubleMatrix(this) - private fun T.wrapVector() = EjmlDoubleVector(this) - - override fun Matrix.unaryMinus(): Matrix = this * elementAlgebra { -one } - - override fun Matrix.dot(other: Matrix): EjmlDoubleMatrix { - val out = DMatrixSparseCSC(1, 1) - CommonOps_DSCC.mult(toEjml().origin, other.toEjml().origin, out) - return out.wrapMatrix() - } - - override fun Matrix.dot(vector: Point): EjmlDoubleVector { - val out = DMatrixSparseCSC(1, 1) - CommonOps_DSCC.mult(toEjml().origin, vector.toEjml().origin, out) - return out.wrapVector() - } - - override operator fun Matrix.minus(other: Matrix): EjmlDoubleMatrix { - val out = DMatrixSparseCSC(1, 1) - - CommonOps_DSCC.add( - elementAlgebra.one, - toEjml().origin, - elementAlgebra { -one }, - other.toEjml().origin, - out, - null, - null, - ) - - return out.wrapMatrix() - } - - override operator fun Matrix.times(value: Double): EjmlDoubleMatrix { - val res = DMatrixSparseCSC(1, 1) - CommonOps_DSCC.scale(value, toEjml().origin, res) - return res.wrapMatrix() - } - - override fun Point.unaryMinus(): EjmlDoubleVector { - val res = DMatrixSparseCSC(1, 1) - CommonOps_DSCC.changeSign(toEjml().origin, res) - return res.wrapVector() - } - - override fun Matrix.plus(other: Matrix): EjmlDoubleMatrix { - val out = DMatrixSparseCSC(1, 1) - - CommonOps_DSCC.add( - elementAlgebra.one, - toEjml().origin, - elementAlgebra.one, - other.toEjml().origin, - out, - null, - null, - ) - - return out.wrapMatrix() - } - - override fun Point.plus(other: Point): EjmlDoubleVector { - val out = DMatrixSparseCSC(1, 1) - - CommonOps_DSCC.add( - elementAlgebra.one, - toEjml().origin, - elementAlgebra.one, - other.toEjml().origin, - out, - null, - null, - ) - - return out.wrapVector() - } - - override fun Point.minus(other: Point): EjmlDoubleVector { - val out = DMatrixSparseCSC(1, 1) - - CommonOps_DSCC.add( - elementAlgebra.one, - toEjml().origin, - elementAlgebra { -one }, - other.toEjml().origin, - out, - null, - null, - ) - - return out.wrapVector() - } - - override fun Double.times(m: Matrix): EjmlDoubleMatrix = m * this - - override fun Point.times(value: Double): EjmlDoubleVector { - val res = DMatrixSparseCSC(1, 1) - CommonOps_DSCC.scale(value, toEjml().origin, res) - return res.wrapVector() - } - - override fun Double.times(v: Point): EjmlDoubleVector = v * this - - @UnstableKMathAPI - override fun computeFeature(structure: Matrix, type: KClass): F? { - structure.getFeature(type)?.let { return it } - val origin = structure.toEjml().origin - - return when (type) { - QRDecompositionFeature::class -> object : QRDecompositionFeature { - private val qr by lazy { - DecompositionFactory_DSCC.qr(FillReducing.NONE).apply { decompose(origin.copy()) } - } - - override val q: Matrix by lazy { - qr.getQ(null, false).wrapMatrix().withFeature(OrthogonalFeature) - } - - override val r: Matrix by lazy { qr.getR(null, false).wrapMatrix().withFeature(UFeature) } - } - - CholeskyDecompositionFeature::class -> object : CholeskyDecompositionFeature { - override val l: Matrix by lazy { - val cholesky = - DecompositionFactory_DSCC.cholesky().apply { decompose(origin.copy()) } - - (cholesky.getT(null) as DMatrix).wrapMatrix().withFeature(LFeature) - } - } - - LUDecompositionFeature::class, DeterminantFeature::class, InverseMatrixFeature::class -> object : - LUDecompositionFeature, DeterminantFeature, InverseMatrixFeature { - private val lu by lazy { - DecompositionFactory_DSCC.lu(FillReducing.NONE).apply { decompose(origin.copy()) } - } - - override val l: Matrix by lazy { - lu.getLower(null).wrapMatrix().withFeature(LFeature) - } - - override val u: Matrix by lazy { - lu.getUpper(null).wrapMatrix().withFeature(UFeature) - } - - override val inverse: Matrix 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(it) - } - } - - /** - * Solves for *x* in the following equation: *x = [a] -1 · [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, b: Matrix): EjmlDoubleMatrix { - 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] -1 · [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, b: Point): EjmlDoubleVector { - 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() { - /** - * The [FloatField] reference. - */ - override val elementAlgebra: FloatField get() = FloatField - - @Suppress("UNCHECKED_CAST") - override fun Matrix.toEjml(): EjmlFloatMatrix = when { - this is EjmlFloatMatrix<*> && origin is FMatrixSparseCSC -> this as EjmlFloatMatrix - else -> buildMatrix(rowNum, colNum) { i, j -> get(i, j) } - } - - @Suppress("UNCHECKED_CAST") - override fun Point.toEjml(): EjmlFloatVector = when { - this is EjmlFloatVector<*> && origin is FMatrixSparseCSC -> this as EjmlFloatVector - else -> EjmlFloatVector(FMatrixSparseCSC(size, 1).also { - (0 until it.numRows).forEach { row -> it[row, 0] = get(row) } - }) - } - - override fun buildMatrix( - rows: Int, - columns: Int, - initializer: FloatField.(i: Int, j: Int) -> Float, - ): EjmlFloatMatrix = FMatrixSparseCSC(rows, columns).also { - (0 until rows).forEach { row -> - (0 until columns).forEach { col -> it[row, col] = elementAlgebra.initializer(row, col) } - } - }.wrapMatrix() - - override fun buildVector( - size: Int, - initializer: FloatField.(Int) -> Float, - ): EjmlFloatVector = EjmlFloatVector(FMatrixSparseCSC(size, 1).also { - (0 until it.numRows).forEach { row -> it[row, 0] = elementAlgebra.initializer(row) } - }) - - private fun T.wrapMatrix() = EjmlFloatMatrix(this) - private fun T.wrapVector() = EjmlFloatVector(this) - - override fun Matrix.unaryMinus(): Matrix = this * elementAlgebra { -one } - - override fun Matrix.dot(other: Matrix): EjmlFloatMatrix { - val out = FMatrixSparseCSC(1, 1) - CommonOps_FSCC.mult(toEjml().origin, other.toEjml().origin, out) - return out.wrapMatrix() - } - - override fun Matrix.dot(vector: Point): EjmlFloatVector { - val out = FMatrixSparseCSC(1, 1) - CommonOps_FSCC.mult(toEjml().origin, vector.toEjml().origin, out) - return out.wrapVector() - } - - override operator fun Matrix.minus(other: Matrix): EjmlFloatMatrix { - val out = FMatrixSparseCSC(1, 1) - - CommonOps_FSCC.add( - elementAlgebra.one, - toEjml().origin, - elementAlgebra { -one }, - other.toEjml().origin, - out, - null, - null, - ) - - return out.wrapMatrix() - } - - override operator fun Matrix.times(value: Float): EjmlFloatMatrix { - val res = FMatrixSparseCSC(1, 1) - CommonOps_FSCC.scale(value, toEjml().origin, res) - return res.wrapMatrix() - } - - override fun Point.unaryMinus(): EjmlFloatVector { - val res = FMatrixSparseCSC(1, 1) - CommonOps_FSCC.changeSign(toEjml().origin, res) - return res.wrapVector() - } - - override fun Matrix.plus(other: Matrix): EjmlFloatMatrix { - val out = FMatrixSparseCSC(1, 1) - - CommonOps_FSCC.add( - elementAlgebra.one, - toEjml().origin, - elementAlgebra.one, - other.toEjml().origin, - out, - null, - null, - ) - - return out.wrapMatrix() - } - - override fun Point.plus(other: Point): EjmlFloatVector { - val out = FMatrixSparseCSC(1, 1) - - CommonOps_FSCC.add( - elementAlgebra.one, - toEjml().origin, - elementAlgebra.one, - other.toEjml().origin, - out, - null, - null, - ) - - return out.wrapVector() - } - - override fun Point.minus(other: Point): EjmlFloatVector { - val out = FMatrixSparseCSC(1, 1) - - CommonOps_FSCC.add( - elementAlgebra.one, - toEjml().origin, - elementAlgebra { -one }, - other.toEjml().origin, - out, - null, - null, - ) - - return out.wrapVector() - } - - override fun Float.times(m: Matrix): EjmlFloatMatrix = m * this - - override fun Point.times(value: Float): EjmlFloatVector { - val res = FMatrixSparseCSC(1, 1) - CommonOps_FSCC.scale(value, toEjml().origin, res) - return res.wrapVector() - } - - override fun Float.times(v: Point): EjmlFloatVector = v * this - - @UnstableKMathAPI - override fun computeFeature(structure: Matrix, type: KClass): F? { - structure.getFeature(type)?.let { return it } - val origin = structure.toEjml().origin - - return when (type) { - QRDecompositionFeature::class -> object : QRDecompositionFeature { - private val qr by lazy { - DecompositionFactory_FSCC.qr(FillReducing.NONE).apply { decompose(origin.copy()) } - } - - override val q: Matrix by lazy { - qr.getQ(null, false).wrapMatrix().withFeature(OrthogonalFeature) - } - - override val r: Matrix by lazy { qr.getR(null, false).wrapMatrix().withFeature(UFeature) } - } - - CholeskyDecompositionFeature::class -> object : CholeskyDecompositionFeature { - override val l: Matrix by lazy { - val cholesky = - DecompositionFactory_FSCC.cholesky().apply { decompose(origin.copy()) } - - (cholesky.getT(null) as FMatrix).wrapMatrix().withFeature(LFeature) - } - } - - LUDecompositionFeature::class, DeterminantFeature::class, InverseMatrixFeature::class -> object : - LUDecompositionFeature, DeterminantFeature, InverseMatrixFeature { - private val lu by lazy { - DecompositionFactory_FSCC.lu(FillReducing.NONE).apply { decompose(origin.copy()) } - } - - override val l: Matrix by lazy { - lu.getLower(null).wrapMatrix().withFeature(LFeature) - } - - override val u: Matrix by lazy { - lu.getUpper(null).wrapMatrix().withFeature(UFeature) - } - - override val inverse: Matrix 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(it) - } - } - - /** - * Solves for *x* in the following equation: *x = [a] -1 · [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, b: Matrix): EjmlFloatMatrix { - 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] -1 · [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, b: Point): EjmlFloatVector { - val res = FMatrixSparseCSC(1, 1) - CommonOps_FSCC.solve(FMatrixSparseCSC(a.toEjml().origin), FMatrixSparseCSC(b.toEjml().origin), res) - return EjmlFloatVector(res) - } -} - diff --git a/kmath-histograms/src/commonMain/kotlin/space/kscience/kmath/histogram/UniformHistogramGroupND.kt b/kmath-histograms/src/commonMain/kotlin/space/kscience/kmath/histogram/UniformHistogramGroupND.kt index 90ec29ce3..eafd55513 100644 --- a/kmath-histograms/src/commonMain/kotlin/space/kscience/kmath/histogram/UniformHistogramGroupND.kt +++ b/kmath-histograms/src/commonMain/kotlin/space/kscience/kmath/histogram/UniformHistogramGroupND.kt @@ -28,7 +28,7 @@ public class UniformHistogramGroupND>( private val lower: Buffer, private val upper: Buffer, private val binNums: IntArray = IntArray(lower.size) { 20 }, - private val bufferFactory: BufferFactory = Buffer.Companion::boxing, + private val bufferFactory: BufferFactory = BufferFactory(Buffer.Companion::boxing), ) : HistogramGroupND { init { @@ -114,7 +114,7 @@ public class UniformHistogramGroupND>( public fun > Histogram.Companion.uniformNDFromRanges( valueAlgebraND: FieldOpsND, vararg ranges: ClosedFloatingPointRange, - bufferFactory: BufferFactory = Buffer.Companion::boxing, + bufferFactory: BufferFactory = BufferFactory(Buffer.Companion::boxing), ): UniformHistogramGroupND = UniformHistogramGroupND( valueAlgebraND, ranges.map(ClosedFloatingPointRange::start).asBuffer(), @@ -140,7 +140,7 @@ public fun Histogram.Companion.uniformDoubleNDFromRanges( public fun > Histogram.Companion.uniformNDFromRanges( valueAlgebraND: FieldOpsND, vararg ranges: Pair, Int>, - bufferFactory: BufferFactory = Buffer.Companion::boxing, + bufferFactory: BufferFactory = BufferFactory(Buffer.Companion::boxing), ): UniformHistogramGroupND = UniformHistogramGroupND( valueAlgebraND, ListBuffer( diff --git a/kmath-multik/src/main/kotlin/space/kscience/kmath/multik/MultikDoubleAlgebra.kt b/kmath-multik/src/main/kotlin/space/kscience/kmath/multik/MultikDoubleAlgebra.kt index 1dc318517..0de2d8349 100644 --- a/kmath-multik/src/main/kotlin/space/kscience/kmath/multik/MultikDoubleAlgebra.kt +++ b/kmath-multik/src/main/kotlin/space/kscience/kmath/multik/MultikDoubleAlgebra.kt @@ -6,6 +6,7 @@ package space.kscience.kmath.multik import org.jetbrains.kotlinx.multik.ndarray.data.DataType +import space.kscience.kmath.misc.PerformancePitfall import space.kscience.kmath.nd.StructureND import space.kscience.kmath.operations.DoubleField import space.kscience.kmath.operations.ExponentialOperations @@ -22,10 +23,13 @@ public object MultikDoubleAlgebra : MultikDivisionTensorAlgebra): MultikTensor = sin(arg) / cos(arg) + @PerformancePitfall override fun asin(arg: StructureND): MultikTensor = arg.map { asin(it) } + @PerformancePitfall override fun acos(arg: StructureND): MultikTensor = arg.map { acos(it) } + @PerformancePitfall override fun atan(arg: StructureND): MultikTensor = arg.map { atan(it) } override fun exp(arg: StructureND): MultikTensor = multikMath.mathEx.exp(arg.asMultik().array).wrap() @@ -42,10 +46,13 @@ public object MultikDoubleAlgebra : MultikDivisionTensorAlgebra): MultikTensor = arg.map { asinh(it) } + @PerformancePitfall override fun acosh(arg: StructureND): MultikTensor = arg.map { acosh(it) } + @PerformancePitfall override fun atanh(arg: StructureND): MultikTensor = arg.map { atanh(it) } } diff --git a/kmath-stat/src/commonMain/kotlin/space/kscience/kmath/stat/Sampler.kt b/kmath-stat/src/commonMain/kotlin/space/kscience/kmath/stat/Sampler.kt index a88f3e437..890318e31 100644 --- a/kmath-stat/src/commonMain/kotlin/space/kscience/kmath/stat/Sampler.kt +++ b/kmath-stat/src/commonMain/kotlin/space/kscience/kmath/stat/Sampler.kt @@ -35,7 +35,7 @@ public fun interface Sampler { public fun Sampler.sampleBuffer( generator: RandomGenerator, size: Int, - bufferFactory: BufferFactory = Buffer.Companion::boxing, + bufferFactory: BufferFactory = BufferFactory(Buffer.Companion::boxing), ): Chain> { require(size > 1) //creating temporary storage once