Finishing fixes

CHANGELOG.md

@@ -20,9 +20,11 @@
 - Kmath-memory is moved on top of core.
 
 ### Deprecated
+- ND4J engine
 
 ### Removed
 - `asPolynomial` function due to scope pollution
+- Codegend for ejml (450 lines of codegen for 1000 lines of code is too much)
 
 ### Fixed
 - Median statistics

benchmarks/src/jsMain/kotlin/space/kscience/kmath/benchmarks/ExpressionsInterpretersBenchmark.kt

@@ -99,7 +99,7 @@ class ExpressionsInterpretersBenchmark {
         private val estree = node.estreeCompileToExpression(Float64Field)
 
         private val raw = Expression<Double> { args ->
-            val x = args[x]!!
+            val x = args.getValue(x)
             x * 2.0 + 2.0 / x - 16.0 / sin(x)
         }
     }

examples/src/main/kotlin/space/kscience/kmath/fit/chiSquared.kt

@@ -13,10 +13,7 @@ import space.kscience.kmath.expressions.autodiff
 import space.kscience.kmath.expressions.symbol
 import space.kscience.kmath.operations.asIterable
 import space.kscience.kmath.operations.toList
-import space.kscience.kmath.optimization.FunctionOptimizationTarget
+import space.kscience.kmath.optimization.*
-import space.kscience.kmath.optimization.optimizeWith
-import space.kscience.kmath.optimization.result
-import space.kscience.kmath.optimization.resultValue
 import space.kscience.kmath.random.RandomGenerator
 import space.kscience.kmath.real.DoubleVector
 import space.kscience.kmath.real.map
@@ -80,8 +77,9 @@ suspend fun main() {
     val result = chi2.optimizeWith(
         CMOptimizer,
         mapOf(a to 1.5, b to 0.9, c to 1.0),
-        FunctionOptimizationTarget.MINIMIZE
+    ){
-    )
+        FunctionOptimizationTarget(OptimizationDirection.MINIMIZE)
+    }
 
     //display a page with plot and numerical results
     val page = Plotly.page {

examples/src/main/kotlin/space/kscience/kmath/fit/qowFit.kt

@@ -7,6 +7,7 @@ package space.kscience.kmath.fit
 
 import kotlinx.html.br
 import kotlinx.html.h3
+import space.kscience.attributes.Attributes
 import space.kscience.kmath.data.XYErrorColumnarData
 import space.kscience.kmath.distributions.NormalDistribution
 import space.kscience.kmath.expressions.Symbol
@@ -64,7 +65,7 @@ suspend fun main() {
         QowOptimizer,
         Double.autodiff,
         mapOf(a to 0.9, b to 1.2, c to 2.0, e to 1.0, d to 1.0, e to 0.0),
-        OptimizationParameters(a, b, c, d)
+        attributes = Attributes(OptimizationParameters, listOf(a, b, c, d))
     ) { arg ->
         //bind variables to autodiff context
         val a by binding

examples/src/main/kotlin/space/kscience/kmath/functions/interpolate.kt

@@ -8,7 +8,6 @@ package space.kscience.kmath.functions
 import space.kscience.kmath.interpolation.SplineInterpolator
 import space.kscience.kmath.interpolation.interpolatePolynomials
 import space.kscience.kmath.operations.Float64Field
-import space.kscience.kmath.structures.Float64Buffer
 import space.kscience.plotly.Plotly
 import space.kscience.plotly.UnstablePlotlyAPI
 import space.kscience.plotly.makeFile
@@ -24,9 +23,7 @@ fun main() {
         x to sin(x)
     }
 
-    val polynomial: PiecewisePolynomial<Double> = SplineInterpolator(
+    val polynomial: PiecewisePolynomial<Double> = SplineInterpolator(Float64Field).interpolatePolynomials(data)
-        Float64Field, ::Float64Buffer
-    ).interpolatePolynomials(data)
 
     val function = polynomial.asFunction(Float64Field, 0.0)
 

examples/src/main/kotlin/space/kscience/kmath/structures/ComplexND.kt

@@ -8,7 +8,7 @@
 package space.kscience.kmath.structures
 
 import space.kscience.kmath.complex.*
-import space.kscience.kmath.linear.transpose
+import space.kscience.kmath.linear.transposed
 import space.kscience.kmath.nd.StructureND
 import space.kscience.kmath.nd.as2D
 import space.kscience.kmath.nd.ndAlgebra
@@ -60,7 +60,7 @@ fun complexExample() {
             val sum = matrix + x + 1.0
 
             //Represent the sum as 2d-structure and transpose
-            sum.as2D().transpose()
+            sum.as2D().transposed()
         }
     }
 }

examples/src/main/kotlin/space/kscience/kmath/structures/StreamDoubleFieldND.kt

@@ -7,6 +7,7 @@ package space.kscience.kmath.structures
 
 import space.kscience.kmath.PerformancePitfall
 import space.kscience.kmath.nd.*
+import space.kscience.kmath.operations.DoubleField
 import space.kscience.kmath.operations.ExtendedField
 import space.kscience.kmath.operations.Float64Field
 import space.kscience.kmath.operations.NumbersAddOps

examples/src/main/kotlin/space/kscience/kmath/structures/StructureWriteBenchmark.kt

@@ -6,20 +6,18 @@
 package space.kscience.kmath.structures
 
 import space.kscience.kmath.nd.BufferND
-import space.kscience.kmath.nd.ShapeND
-import space.kscience.kmath.nd.StructureND
 import space.kscience.kmath.operations.mapToBuffer
 import kotlin.system.measureTimeMillis
 
 private inline fun <T, reified R : Any> BufferND<T>.mapToBufferND(
-    bufferFactory: BufferFactory<R> = BufferFactory.auto(),
+    bufferFactory: BufferFactory<R> = BufferFactory(),
     crossinline block: (T) -> R,
 ): BufferND<R> = BufferND(indices, buffer.mapToBuffer(bufferFactory, block))
 
 @Suppress("UNUSED_VARIABLE")
 fun main() {
     val n = 6000
-    val structure = StructureND.buffered(ShapeND(n, n), Buffer.Companion::auto) { 1.0 }
+    val structure = BufferND(n, n) { 1.0 }
     structure.mapToBufferND { it + 1 } // warm-up
     val time1 = measureTimeMillis { val res = structure.mapToBufferND { it + 1 } }
     println("Structure mapping finished in $time1 millis")

examples/src/main/kotlin/space/kscience/kmath/structures/buffers.kt

@@ -13,7 +13,7 @@ import space.kscience.kmath.operations.withSize
 inline fun <reified R : Any> MutableBuffer.Companion.same(
     n: Int,
     value: R
-): MutableBuffer<R> = auto(n) { value }
+): MutableBuffer<R> = MutableBuffer(n) { value }
 
 
 fun main() {

kmath-commons/src/jvmMain/kotlin/space/kscience/kmath/commons/linear/CMMatrix.kt

@@ -7,16 +7,20 @@ package space.kscience.kmath.commons.linear
 
 import org.apache.commons.math3.linear.*
 import org.apache.commons.math3.linear.LUDecomposition
+import org.apache.commons.math3.linear.SingularValueDecomposition
 import space.kscience.attributes.SafeType
 import space.kscience.kmath.UnstableKMathAPI
 import space.kscience.kmath.linear.*
+import space.kscience.kmath.linear.CholeskyDecomposition
+import space.kscience.kmath.linear.QRDecomposition
 import space.kscience.kmath.nd.Structure2D
 import space.kscience.kmath.nd.StructureAttribute
 import space.kscience.kmath.operations.DoubleField
 import space.kscience.kmath.operations.Float64Field
 import space.kscience.kmath.structures.Buffer
-import space.kscience.kmath.structures.Float64Buffer
+import space.kscience.kmath.structures.Float64
-import kotlin.reflect.cast
+import space.kscience.kmath.structures.IntBuffer
+import space.kscience.kmath.structures.asBuffer
 
 public class CMMatrix(public val origin: RealMatrix) : Matrix<Double> {
     override val type: SafeType<Double> get() = DoubleField.type
@@ -109,45 +113,44 @@ public object CMLinearSpace : LinearSpace<Double, Float64Field> {
 
         val origin = structure.toCM().origin
 
-        return when (attribute) {
+        val raw: Any? = when (attribute) {
-            IsDiagonal -> if (origin is DiagonalMatrix) IsDiagonal else null
+            IsDiagonal -> if (origin is DiagonalMatrix) Unit else null
             Determinant -> LUDecomposition(origin).determinant
-            LUP -> GenericLupDecomposition {
+
-                private val lup by lazy { LUDecomposition(origin) }
+            LUP -> object : LupDecomposition<Float64> {
-                override val determinant: Double by lazy { lup.determinant }
+                val lup by lazy { LUDecomposition(origin) }
-                override val l: Matrix<Double> by lazy<Matrix<Double>> { CMMatrix(lup.l).withAttribute(LowerTriangular) }
+                override val pivot: IntBuffer get() = lup.pivot.asBuffer()
-                override val u: Matrix<Double> by lazy<Matrix<Double>> { CMMatrix(lup.u).withAttribute(UpperTriangular) }
+                override val l: Matrix<Float64> get() = lup.l.wrap()
-                override val p: Matrix<Double> by lazy { CMMatrix(lup.p) }
+                override val u: Matrix<Float64> get() = lup.u.wrap()
             }
 
-            CholeskyDecompositionAttribute -> object : CholeskyDecompositionAttribute<Double> {
+            Cholesky -> object : CholeskyDecomposition<Float64> {
-                override val l: Matrix<Double> by lazy<Matrix<Double>> {
+                val cmCholesky by lazy { org.apache.commons.math3.linear.CholeskyDecomposition(origin) }
-                    val cholesky = CholeskyDecomposition(origin)
+                override val l: Matrix<Double> get() = cmCholesky.l.wrap()
-                    CMMatrix(cholesky.l).withAttribute(LowerTriangular)
-                }
             }
 
-            QRDecompositionAttribute -> object : QRDecompositionAttribute<Double> {
+            QR -> object : QRDecomposition<Float64> {
-                private val qr by lazy { QRDecomposition(origin) }
+                val cmQr by lazy { org.apache.commons.math3.linear.QRDecomposition(origin) }
-                override val q: Matrix<Double> by lazy<Matrix<Double>> {
+                override val q: Matrix<Float64> get() = cmQr.q.wrap().withAttribute(OrthogonalAttribute)
-                    CMMatrix(qr.q).withAttribute(
+                override val r: Matrix<Float64> get() = cmQr.r.wrap().withAttribute(UpperTriangular)
-                        OrthogonalAttribute
-                    )
-                }
-                override val r: Matrix<Double> by lazy<Matrix<Double>> { CMMatrix(qr.r).withAttribute(UpperTriangular) }
             }
 
-            SVDAttribute -> object : SVDAttribute<Double> {
+            SVD -> object : space.kscience.kmath.linear.SingularValueDecomposition<Float64> {
-                private val sv by lazy { SingularValueDecomposition(origin) }
+                val cmSvd by lazy { SingularValueDecomposition(origin) }
-                override val u: Matrix<Double> by lazy { CMMatrix(sv.u) }
+
-                override val s: Matrix<Double> by lazy { CMMatrix(sv.s) }
+                override val u: Matrix<Float64> get() = cmSvd.u.wrap()
-                override val v: Matrix<Double> by lazy { CMMatrix(sv.v) }
+                override val s: Matrix<Float64> get() = cmSvd.s.wrap()
-                override val singularValues: Point<Double> by lazy { Float64Buffer(sv.singularValues) }
+                override val v: Matrix<Float64> get() = cmSvd.v.wrap()
+                override val singularValues: Point<Float64> get() = cmSvd.singularValues.asBuffer()
+
             }
 
             else -> null
-        }?.let(type::cast)
+        }
+        @Suppress("UNCHECKED_CAST")
+        return raw as V?
     }
+
 }
 
 public operator fun CMMatrix.plus(other: CMMatrix): CMMatrix = CMMatrix(origin.add(other.origin))

kmath-commons/src/jvmMain/kotlin/space/kscience/kmath/commons/optimization/CMOptimizer.kt

@@ -3,6 +3,7 @@
  * Use of this source code is governed by the Apache 2.0 license that can be found in the license/LICENSE.txt file.
  */
 @file:OptIn(UnstableKMathAPI::class)
+
 package space.kscience.kmath.commons.optimization
 
 import org.apache.commons.math3.optim.*
@@ -20,7 +21,6 @@ import space.kscience.kmath.expressions.Symbol
 import space.kscience.kmath.expressions.SymbolIndexer
 import space.kscience.kmath.expressions.derivative
 import space.kscience.kmath.expressions.withSymbols
-import space.kscience.kmath.misc.log
 import space.kscience.kmath.optimization.*
 import kotlin.collections.set
 import kotlin.reflect.KClass
@@ -28,7 +28,7 @@ import kotlin.reflect.KClass
 public operator fun PointValuePair.component1(): DoubleArray = point
 public operator fun PointValuePair.component2(): Double = value
 
-public object CMOptimizerEngine: OptimizationAttribute<() -> MultivariateOptimizer>
+public object CMOptimizerEngine : OptimizationAttribute<() -> MultivariateOptimizer>
 
 /**
  * Specify a Commons-maths optimization engine
@@ -37,7 +37,7 @@ public fun AttributesBuilder<FunctionOptimization<Double>>.cmEngine(optimizerBui
     set(CMOptimizerEngine, optimizerBuilder)
 }
 
-public object CMOptimizerData: SetAttribute<SymbolIndexer.() -> OptimizationData>
+public object CMOptimizerData : SetAttribute<SymbolIndexer.() -> OptimizationData>
 
 /**
  * Specify Commons-maths optimization data.
@@ -118,21 +118,24 @@ public object CMOptimizer : Optimizer<Double, FunctionOptimization<Double>> {
 
             val logger = problem.attributes[OptimizationLog]
 
-            for (feature in problem.attributes) {
+            problem.attributes[CMOptimizerData]?.let { builders: Set<SymbolIndexer.() -> OptimizationData> ->
-                when (feature) {
+                builders.forEach { dataBuilder ->
-                    is CMOptimizerData -> feature.data.forEach { dataBuilder ->
+                    addOptimizationData(dataBuilder())
-                        addOptimizationData(dataBuilder())
+                }
-                    }
+            }
-                    is FunctionOptimizationTarget -> when (feature) {
+
-                        FunctionOptimizationTarget.MAXIMIZE -> addOptimizationData(GoalType.MAXIMIZE)
+            problem.attributes[FunctionOptimizationTarget]?.let { direction: OptimizationDirection ->
-                        FunctionOptimizationTarget.MINIMIZE -> addOptimizationData(GoalType.MINIMIZE)
+                when (direction) {
-                    }
+                    OptimizationDirection.MAXIMIZE -> addOptimizationData(GoalType.MAXIMIZE)
-                    else -> logger?.log { "The feature $feature is unused in optimization" }
+                    OptimizationDirection.MINIMIZE -> addOptimizationData(GoalType.MINIMIZE)
                 }
             }
 
             val (point, value) = cmOptimizer.optimize(*optimizationData.values.toTypedArray())
-            return problem.withAttributes(OptimizationResult(point.toMap()), OptimizationValue(value))
+            return problem.withAttributes {
+                result(point.toMap())
+                value(value)
+            }
         }
     }
 }

kmath-commons/src/jvmTest/kotlin/space/kscience/kmath/commons/integration/IntegrationTest.kt

@@ -7,6 +7,8 @@ package space.kscience.kmath.commons.integration
 
 import org.junit.jupiter.api.Test
 import space.kscience.kmath.UnstableKMathAPI
+import space.kscience.kmath.integration.IntegrandAbsoluteAccuracy
+import space.kscience.kmath.integration.IntegrandRelativeAccuracy
 import space.kscience.kmath.integration.integrate
 import space.kscience.kmath.integration.value
 import space.kscience.kmath.operations.Float64Field.sin
@@ -27,8 +29,8 @@ internal class IntegrationTest {
     @Test
     fun customSimpson() {
         val res = CMIntegrator.simpson().integrate(0.0..PI, {
-            targetRelativeAccuracy = 1e-4
+            IntegrandRelativeAccuracy(1e-4)
-            targetAbsoluteAccuracy = 1e-4
+            IntegrandAbsoluteAccuracy(1e-4)
         }, function).value
         assertTrue { abs(res - 2) < 1e-3 }
         assertTrue { abs(res - 2) > 1e-12 }

kmath-commons/src/jvmTest/kotlin/space/kscience/kmath/commons/optimization/OptimizeTest.kt

@@ -73,8 +73,9 @@ internal class OptimizeTest {
         val result: FunctionOptimization<Double> = chi2.optimizeWith(
             CMOptimizer,
             mapOf(a to 1.5, b to 0.9, c to 1.0),
-            FunctionOptimizationTarget.MINIMIZE
+        ){
-        )
+            FunctionOptimizationTarget(OptimizationDirection.MINIMIZE)
+        }
         println(result)
         println("Chi2/dof = ${result.resultValue / (x.size - 3)}")
     }

kmath-core/src/commonMain/kotlin/space/kscience/kmath/expressions/Expression.kt

@@ -7,6 +7,7 @@ package space.kscience.kmath.expressions
 
 import space.kscience.attributes.SafeType
 import space.kscience.attributes.WithType
+import space.kscience.attributes.safeTypeOf
 import space.kscience.kmath.UnstableKMathAPI
 import space.kscience.kmath.operations.Algebra
 import space.kscience.kmath.operations.DoubleField
@@ -30,12 +31,18 @@ public interface Expression<T> : WithType<T> {
     public operator fun invoke(arguments: Map<Symbol, T>): T
 }
 
+/**
+ * Create an expression from a functional block.
+ */
 public fun <T> Expression(type: SafeType<T>, block: (Map<Symbol, T>) -> T): Expression<T> = object : Expression<T> {
     override fun invoke(arguments: Map<Symbol, T>): T = block(arguments)
 
     override val type: SafeType<T> = type
 }
 
+public inline fun <reified T> Expression(noinline block: (Map<Symbol, T>) -> T): Expression<T> =
+    Expression(safeTypeOf<T>(), block)
+
 /**
  * Specialization of [Expression] for [Double] allowing better performance because of using array.
  */

kmath-core/src/commonMain/kotlin/space/kscience/kmath/linear/LupDecomposition.kt

@@ -15,14 +15,19 @@ import space.kscience.kmath.operations.*
 import space.kscience.kmath.structures.*
 
 public interface LupDecomposition<T> {
-    public val linearSpace: LinearSpace<T, Field<T>>
-    public val elementAlgebra: Field<T> get() = linearSpace.elementAlgebra
-
     public val pivot: IntBuffer
     public val l: Matrix<T>
     public val u: Matrix<T>
 }
 
+/**
+ * Create a pivot matrix from pivot vector using provided [LinearSpace]
+ */
+public fun <T> LupDecomposition<T>.pivotMatrix(linearSpace: LinearSpace<T, Ring<T>>): Matrix<T> =
+    VirtualMatrix(linearSpace.type, l.rowNum, l.colNum) { row, column ->
+        if (column == pivot[row]) linearSpace.elementAlgebra.one else linearSpace.elementAlgebra.zero
+    }
+
 /**
  * Matrices with this feature support LU factorization with partial pivoting: *[p] &middot; a = [l] &middot; [u]* where
  * *a* is the owning matrix.
@@ -31,12 +36,14 @@ public interface LupDecomposition<T> {
  * @param lu combined L and U matrix
  */
 public class GenericLupDecomposition<T>(
-    override val linearSpace: LinearSpace<T, Field<T>>,
+    public val linearSpace: LinearSpace<T, Field<T>>,
     private val lu: Matrix<T>,
     override val pivot: IntBuffer,
     private val even: Boolean,
 ) : LupDecomposition<T> {
 
+    private val elementAlgebra get() = linearSpace.elementAlgebra
+
     override val l: Matrix<T>
         get() = VirtualMatrix(lu.type, lu.rowNum, lu.colNum, attributes = Attributes(LowerTriangular)) { i, j ->
             when {
@@ -51,11 +58,6 @@ public class GenericLupDecomposition<T>(
             if (j >= i) lu[i, j] else elementAlgebra.zero
         }
 
-    public val pivotMatrix: Matrix<T>
-        get() = VirtualMatrix(linearSpace.type, l.rowNum, l.colNum) { row, column ->
-            if (column == pivot[row]) elementAlgebra.one else elementAlgebra.zero
-        }
-
     public val determinant: T by lazy {
         elementAlgebra { (0 until l.shape[0]).fold(if (even) one else -one) { value, i -> value * lu[i, i] } }
     }
@@ -79,7 +81,7 @@ internal fun <T : Comparable<T>> LinearSpace<T, Ring<T>>.abs(value: T): T =
 public fun <T : Comparable<T>> LinearSpace<T, Field<T>>.lup(
     matrix: Matrix<T>,
     checkSingular: (T) -> Boolean,
-): LupDecomposition<T> = elementAlgebra {
+): GenericLupDecomposition<T> = elementAlgebra {
     require(matrix.rowNum == matrix.colNum) { "LU decomposition supports only square matrices" }
     val m = matrix.colNum
     val pivot = IntArray(matrix.rowNum)
@@ -156,7 +158,7 @@ public fun <T : Comparable<T>> LinearSpace<T, Field<T>>.lup(
 public fun LinearSpace<Double, Float64Field>.lup(
     matrix: Matrix<Double>,
     singularityThreshold: Double = 1e-11,
-): LupDecomposition<Double> = lup(matrix) { it < singularityThreshold }
+): GenericLupDecomposition<Double> = lup(matrix) { it < singularityThreshold }
 
 internal fun <T> LinearSpace<T, Field<T>>.solve(
     lup: LupDecomposition<T>,

kmath-core/src/commonMain/kotlin/space/kscience/kmath/linear/Transposed.kt

@@ -25,5 +25,4 @@ public class TransposedMatrix<T>(public val origin: Matrix<T>) : Matrix<T> {
 /**
  * Create a virtual transposed matrix without copying anything. `A.transpose().transpose() === A`
  */
-public val <T> Matrix<T>.transposed: Matrix<T>
+public fun <T> Matrix<T>.transposed(): Matrix<T> = (this as? TransposedMatrix<T>)?.origin ?: TransposedMatrix(this)
-    get() = (this as? TransposedMatrix<T>)?.origin ?: TransposedMatrix(this)

kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/Float64FieldND.kt

@@ -5,9 +5,11 @@
 
 package space.kscience.kmath.nd
 
+import space.kscience.attributes.SafeType
 import space.kscience.kmath.PerformancePitfall
 import space.kscience.kmath.UnstableKMathAPI
 import space.kscience.kmath.operations.*
+import space.kscience.kmath.structures.Float64
 import space.kscience.kmath.structures.Float64Buffer
 import kotlin.contracts.InvocationKind
 import kotlin.contracts.contract
@@ -20,7 +22,10 @@ import kotlin.math.pow as kpow
 public class Float64BufferND(
     indexes: ShapeIndexer,
     override val buffer: Float64Buffer,
-) : MutableBufferND<Double>(indexes, buffer), MutableStructureNDOfDouble{
+) : MutableBufferND<Double>(indexes, buffer), MutableStructureNDOfDouble {
+
+    override val type: SafeType<Float64> get() = Float64Field.type
+
     override fun getDouble(index: IntArray): Double = buffer[indices.offset(index)]
 
     override fun setDouble(index: IntArray, value: Double) {

kmath-core/src/commonMain/kotlin/space/kscience/kmath/nd/primitiveStructureND.kt

@@ -5,9 +5,15 @@
 
 package space.kscience.kmath.nd
 
+import space.kscience.attributes.SafeType
 import space.kscience.kmath.PerformancePitfall
+import space.kscience.kmath.operations.Float64Field
+import space.kscience.kmath.operations.Int32Field
+import space.kscience.kmath.structures.Float64
 
 public interface StructureNDOfDouble : StructureND<Double> {
+    override val type: SafeType<Float64> get() = Float64Field.type
+
     /**
      * Guaranteed non-blocking access to content
      */
@@ -22,6 +28,7 @@ public fun StructureND<Double>.getDouble(index: IntArray): Double =
     if (this is StructureNDOfDouble) getDouble(index) else get(index)
 
 public interface MutableStructureNDOfDouble : StructureNDOfDouble, MutableStructureND<Double> {
+
     /**
      * Guaranteed non-blocking access to content
      */
@@ -34,6 +41,9 @@ public fun MutableStructureND<Double>.getDouble(index: IntArray): Double =
 
 
 public interface StructureNDOfInt : StructureND<Int> {
+
+    override val type: SafeType<Int> get() = Int32Field.type
+
     /**
      * Guaranteed non-blocking access to content
      */

kmath-core/src/commonTest/kotlin/space/kscience/kmath/linear/DoubleLUSolverTest.kt

@@ -40,7 +40,7 @@ class DoubleLUSolverTest {
         //Check determinant
         assertEquals(7.0, lup.determinant)
 
-        assertMatrixEquals(lup.pivotMatrix dot matrix, lup.l dot lup.u)
+        assertMatrixEquals(lup.pivotMatrix(this) dot matrix, lup.l dot lup.u)
     }
 
     @Test

kmath-core/src/commonTest/kotlin/space/kscience/kmath/linear/MatrixTest.kt

@@ -21,7 +21,7 @@ class MatrixTest {
     @Test
     fun testTranspose() = Double.algebra.linearSpace.run {
         val matrix = one(3, 3)
-        val transposed = matrix.transposed
+        val transposed = matrix.transposed()
         assertTrue { StructureND.contentEquals(matrix, transposed) }
     }
 

kmath-coroutines/build.gradle.kts

@@ -6,6 +6,7 @@ kscience {
     jvm()
     js()
     native()
+    wasm()
 
     dependencies {
         api(project(":kmath-core"))

kmath-coroutines/src/jvmTest/kotlin/space/kscience/kmath/streaming/RingBufferTest.kt

@@ -7,6 +7,7 @@ package space.kscience.kmath.streaming
 
 import kotlinx.coroutines.flow.*
 import kotlinx.coroutines.runBlocking
+import space.kscience.kmath.operations.Int32Ring
 import space.kscience.kmath.operations.asSequence
 import kotlin.test.Test
 import kotlin.test.assertEquals
@@ -14,7 +15,7 @@ import kotlin.test.assertEquals
 internal class RingBufferTest {
     @Test
     fun push() {
-        val buffer = RingBuffer.build(20, Double.NaN)
+        val buffer = RingBuffer(20, Double.NaN)
         runBlocking {
             for (i in 1..30) {
                 buffer.push(i.toDouble())
@@ -30,7 +31,7 @@ internal class RingBufferTest {
             while (true) emit(i++)
         }
 
-        val windowed = flow.windowed(10)
+        val windowed = flow.windowed(10, Int32Ring)
 
         runBlocking {
             @Suppress("UNUSED_VARIABLE") val first = windowed.take(1).single()

kmath-dimensions/build.gradle.kts

@@ -6,6 +6,7 @@ kscience{
     jvm()
     js()
     native()
+    wasm()
 
     dependencies{
         api(projects.kmathCore)

kmath-dimensions/src/commonMain/kotlin/space/kscience/kmath/dimensions/Wrappers.kt

@@ -153,7 +153,7 @@ public value class DMatrixContext<T : Any, out A : Ring<T>>(public val context:
         context.run { this@unaryMinus.unaryMinus() }.coerce()
 
     public inline fun <reified R : Dimension, reified C : Dimension> DMatrix<T, C, R>.transposed(): DMatrix<T, R, C> =
-        context.run { (this@transposed as Matrix<T>).transposed }.coerce()
+        context.run { (this@transposed as Matrix<T>).transposed() }.coerce()
 
     public companion object {
         public val real: DMatrixContext<Double, Float64Field> = DMatrixContext(Double.algebra.linearSpace)

kmath-dimensions/src/wasmJsMain/kotlin/space/kscience/kmath/dimensions/Dimension.wasmJs.kt

@@ -0,0 +1,23 @@
+/*
+ * Copyright 2018-2024 KMath contributors.
+ * Use of this source code is governed by the Apache 2.0 license that can be found in the license/LICENSE.txt file.
+ */
+
+package space.kscience.kmath.dimensions
+
+import kotlin.reflect.KClass
+
+private val dimensionMap: MutableMap<Int, Dimension> = hashMapOf(1 to D1, 2 to D2, 3 to D3)
+
+@Suppress("UNCHECKED_CAST")
+public actual fun <D : Dimension> Dimension.Companion.resolve(type: KClass<D>): D = dimensionMap
+    .entries
+    .map(MutableMap.MutableEntry<Int, Dimension>::value)
+    .find { it::class == type } as? D
+    ?: error("Can't resolve dimension $type")
+
+public actual fun Dimension.Companion.of(dim: Int): Dimension = dimensionMap.getOrPut(dim) {
+    object : Dimension {
+        override val dim: Int get() = dim
+    }
+}

kmath-ejml/build.gradle.kts

@@ -1,15 +1,15 @@
-import space.kscience.kmath.ejml.codegen.ejmlCodegen
-
 plugins {
     id("space.kscience.gradle.jvm")
 }
 
+val ejmlVerision = "0.43.1"
+
 dependencies {
-    api("org.ejml:ejml-ddense:0.41")
+    api("org.ejml:ejml-ddense:$ejmlVerision")
-    api("org.ejml:ejml-fdense:0.41")
+    api("org.ejml:ejml-fdense:$ejmlVerision")
-    api("org.ejml:ejml-dsparse:0.41")
+    api("org.ejml:ejml-dsparse:$ejmlVerision")
-    api("org.ejml:ejml-fsparse:0.41")
+    api("org.ejml:ejml-fsparse:$ejmlVerision")
-    api(project(":kmath-core"))
+    api(projects.kmathCore)
 }
 
 readme {
@@ -32,10 +32,10 @@ readme {
     ) { "LinearSpace implementations." }
 }
 
-kotlin.sourceSets.main {
+//kotlin.sourceSets.main {
-    val codegen by tasks.creating {
+//    val codegen by tasks.creating {
-        ejmlCodegen(kotlin.srcDirs.first().absolutePath + "/space/kscience/kmath/ejml/_generated.kt")
+//        ejmlCodegen(kotlin.srcDirs.first().absolutePath + "/space/kscience/kmath/ejml/_generated.kt")
-    }
+//    }
-
+//
-    kotlin.srcDirs(files().builtBy(codegen))
+//    kotlin.srcDirs(files().builtBy(codegen))
-}
+//}

kmath-ejml/src/main/kotlin/space/kscience/kmath/ejml/implementations.kt

@@ -17,21 +17,18 @@ 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.attributes.SafeType
 import space.kscience.attributes.safeTypeOf
-import space.kscience.kmath.UnstableKMathAPI
 import space.kscience.kmath.linear.*
 import space.kscience.kmath.linear.Matrix
-import space.kscience.kmath.nd.StructureFeature
+import space.kscience.kmath.nd.Structure2D
+import space.kscience.kmath.nd.StructureAttribute
 import space.kscience.kmath.operations.Float32Field
 import space.kscience.kmath.operations.Float64Field
 import space.kscience.kmath.operations.invoke
-import space.kscience.kmath.structures.DoubleBuffer
+import space.kscience.kmath.structures.Float32
-import space.kscience.kmath.structures.FloatBuffer
+import space.kscience.kmath.structures.IntBuffer
-import kotlin.reflect.KClass
+import space.kscience.kmath.structures.asBuffer
-import kotlin.reflect.cast
 
 /**
  * [EjmlVector] specialization for [Double].
@@ -78,7 +75,6 @@ public class EjmlFloatMatrix<out M : FMatrix>(override val origin: M) : EjmlMatr
 }
 
 
-
 /**
  * [EjmlLinearSpace] implementation based on [CommonOps_DDRM], [DecompositionFactory_DDRM] operations and
  * [DMatrixRMaj] matrices.
@@ -217,77 +213,65 @@ public object EjmlLinearSpaceDDRM : EjmlLinearSpace<Double, Float64Field, DMatri
 
     override fun Double.times(v: Point<Double>): EjmlDoubleVector<DMatrixRMaj> = v * this
 
-    @UnstableKMathAPI
+    override fun <V, A : StructureAttribute<V>> computeAttribute(structure: Structure2D<Double>, attribute: A): V? {
-    override fun <F : StructureFeature> computeFeature(structure: Matrix<Double>, type: KClass<out F>): F? {
-        structure.getFeature(type)?.let { return it }
         val origin = structure.toEjml().origin
 
-        return when (type) {
+        val raw: Any? = when (attribute) {
-            InverseMatrixFeature::class -> object : InverseMatrixFeature<Double> {
+            Inverted -> {
-                override val inverse: Matrix<Double> by lazy {
+                val res = origin.copy()
-                    val res = origin.copy()
+                CommonOps_DDRM.invert(res)
-                    CommonOps_DDRM.invert(res)
+                res.wrapMatrix()
-                    res.wrapMatrix()
-                }
             }
 
-            DeterminantFeature::class -> object : DeterminantFeature<Double> {
+            Determinant -> CommonOps_DDRM.det(origin)
-                override val determinant: Double by lazy { CommonOps_DDRM.det(origin) }
+            SVD -> object : SingularValueDecomposition<Double> {
-            }
+                val ejmlSvd by lazy {
-
+                    DecompositionFactory_DDRM
-            SingularValueDecompositionFeature::class -> object : SingularValueDecompositionFeature<Double> {
+                        .svd(origin.numRows, origin.numCols, true, true, false)
-                private val svd by lazy {
-                    DecompositionFactory_DDRM.svd(origin.numRows, origin.numCols, true, true, false)
                         .apply { decompose(origin.copy()) }
                 }
+                override val u: Matrix<Double> get() = ejmlSvd.getU(null, false).wrapMatrix()
+
+                override val s: Matrix<Double> get() = ejmlSvd.getW(null).wrapMatrix()
+                override val v: Matrix<Double> get() = ejmlSvd.getV(null, false).wrapMatrix()
+                override val singularValues: Point<Double> get() = ejmlSvd.singularValues.asBuffer()
 
-                override val u: Matrix<Double> by lazy { svd.getU(null, false).wrapMatrix() }
-                override val s: Matrix<Double> by lazy { svd.getW(null).wrapMatrix() }
-                override val v: Matrix<Double> by lazy { svd.getV(null, false).wrapMatrix() }
-                override val singularValues: Point<Double> by lazy { DoubleBuffer(svd.singularValues) }
             }
 
-            QRDecompositionFeature::class -> object : QRDecompositionFeature<Double> {
+            QR -> object : QRDecomposition<Double> {
-                private val qr by lazy {
+                val ejmlQr by lazy { DecompositionFactory_DDRM.qr().apply { decompose(origin.copy()) } }
-                    DecompositionFactory_DDRM.qr().apply { decompose(origin.copy()) }
+                override val q: Matrix<Double> get() = ejmlQr.getQ(null, false).wrapMatrix()
-                }
+                override val r: Matrix<Double> get() = ejmlQr.getR(null, false).wrapMatrix()
-
-                override val q: Matrix<Double> by lazy {
-                    qr.getQ(null, false).wrapMatrix().withFeature(OrthogonalFeature)
-                }
-
-                override val r: Matrix<Double> by lazy { qr.getR(null, false).wrapMatrix().withFeature(UFeature) }
             }
 
-            CholeskyDecompositionFeature::class -> object : CholeskyDecompositionFeature<Double> {
+            Cholesky -> object : CholeskyDecomposition<Double> {
                 override val l: Matrix<Double> by lazy {
                     val cholesky =
                         DecompositionFactory_DDRM.chol(structure.rowNum, true).apply { decompose(origin.copy()) }
 
-                    cholesky.getT(null).wrapMatrix().withFeature(LFeature)
+                    cholesky.getT(null).wrapMatrix().withAttribute(LowerTriangular)
                 }
             }
 
-            LupDecompositionFeature::class -> object : LupDecompositionFeature<Double> {
+            LUP -> object : LupDecomposition<Double> {
                 private val lup by lazy {
                     DecompositionFactory_DDRM.lu(origin.numRows, origin.numCols).apply { decompose(origin.copy()) }
                 }
 
-                override val l: Matrix<Double> by lazy {
+                override val l: Matrix<Double>
-                    lup.getLower(null).wrapMatrix().withFeature(LFeature)
+                    get() = lup.getLower(null).wrapMatrix().withAttribute(LowerTriangular)
-                }
 
-                override val u: Matrix<Double> by lazy {
-                    lup.getUpper(null).wrapMatrix().withFeature(UFeature)
-                }
 
-                override val p: Matrix<Double> by lazy { lup.getRowPivot(null).wrapMatrix() }
+                override val u: Matrix<Double>
+                    get() = lup.getUpper(null).wrapMatrix().withAttribute(UpperTriangular)
+                override val pivot: IntBuffer get() = lup.getRowPivotV(null).asBuffer()
             }
 
             else -> null
-        }?.let{
-            type.cast(it)
         }
+
+        @Suppress("UNCHECKED_CAST")
+        return raw as V?
     }
 
     /**
@@ -318,7 +302,6 @@ public object EjmlLinearSpaceDDRM : EjmlLinearSpace<Double, Float64Field, DMatri
 }
 
 
-
 /**
  * [EjmlLinearSpace] implementation based on [CommonOps_FDRM], [DecompositionFactory_FDRM] operations and
  * [FMatrixRMaj] matrices.
@@ -457,77 +440,65 @@ public object EjmlLinearSpaceFDRM : EjmlLinearSpace<Float, Float32Field, FMatrix
 
     override fun Float.times(v: Point<Float>): EjmlFloatVector<FMatrixRMaj> = v * this
 
-    @UnstableKMathAPI
+    override fun <V, A : StructureAttribute<V>> computeAttribute(structure: Structure2D<Float32>, attribute: A): V? {
-    override fun <F : StructureFeature> computeFeature(structure: Matrix<Float>, type: KClass<out F>): F? {
-        structure.getFeature(type)?.let { return it }
         val origin = structure.toEjml().origin
 
-        return when (type) {
+        val raw: Any? = when (attribute) {
-            InverseMatrixFeature::class -> object : InverseMatrixFeature<Float> {
+            Inverted -> {
-                override val inverse: Matrix<Float> by lazy {
+                val res = origin.copy()
-                    val res = origin.copy()
+                CommonOps_FDRM.invert(res)
-                    CommonOps_FDRM.invert(res)
+                res.wrapMatrix()
-                    res.wrapMatrix()
-                }
             }
 
-            DeterminantFeature::class -> object : DeterminantFeature<Float> {
+            Determinant -> CommonOps_FDRM.det(origin)
-                override val determinant: Float by lazy { CommonOps_FDRM.det(origin) }
+            SVD -> object : SingularValueDecomposition<Float32> {
-            }
+                val ejmlSvd by lazy {
-
+                    DecompositionFactory_FDRM
-            SingularValueDecompositionFeature::class -> object : SingularValueDecompositionFeature<Float> {
+                        .svd(origin.numRows, origin.numCols, true, true, false)
-                private val svd by lazy {
-                    DecompositionFactory_FDRM.svd(origin.numRows, origin.numCols, true, true, false)
                         .apply { decompose(origin.copy()) }
                 }
+                override val u: Matrix<Float32> get() = ejmlSvd.getU(null, false).wrapMatrix()
+
+                override val s: Matrix<Float32> get() = ejmlSvd.getW(null).wrapMatrix()
+                override val v: Matrix<Float32> get() = ejmlSvd.getV(null, false).wrapMatrix()
+                override val singularValues: Point<Float32> get() = ejmlSvd.singularValues.asBuffer()
 
-                override val u: Matrix<Float> by lazy { svd.getU(null, false).wrapMatrix() }
-                override val s: Matrix<Float> by lazy { svd.getW(null).wrapMatrix() }
-                override val v: Matrix<Float> by lazy { svd.getV(null, false).wrapMatrix() }
-                override val singularValues: Point<Float> by lazy { FloatBuffer(svd.singularValues) }
             }
 
-            QRDecompositionFeature::class -> object : QRDecompositionFeature<Float> {
+            QR -> object : QRDecomposition<Float32> {
-                private val qr by lazy {
+                val ejmlQr by lazy { DecompositionFactory_FDRM.qr().apply { decompose(origin.copy()) } }
-                    DecompositionFactory_FDRM.qr().apply { decompose(origin.copy()) }
+                override val q: Matrix<Float32> get() = ejmlQr.getQ(null, false).wrapMatrix()
-                }
+                override val r: Matrix<Float32> get() = ejmlQr.getR(null, false).wrapMatrix()
-
-                override val q: Matrix<Float> by lazy {
-                    qr.getQ(null, false).wrapMatrix().withFeature(OrthogonalFeature)
-                }
-
-                override val r: Matrix<Float> by lazy { qr.getR(null, false).wrapMatrix().withFeature(UFeature) }
             }
 
-            CholeskyDecompositionFeature::class -> object : CholeskyDecompositionFeature<Float> {
+            Cholesky -> object : CholeskyDecomposition<Float32> {
-                override val l: Matrix<Float> by lazy {
+                override val l: Matrix<Float32> by lazy {
                     val cholesky =
                         DecompositionFactory_FDRM.chol(structure.rowNum, true).apply { decompose(origin.copy()) }
 
-                    cholesky.getT(null).wrapMatrix().withFeature(LFeature)
+                    cholesky.getT(null).wrapMatrix().withAttribute(LowerTriangular)
                 }
             }
 
-            LupDecompositionFeature::class -> object : LupDecompositionFeature<Float> {
+            LUP -> object : LupDecomposition<Float32> {
                 private val lup by lazy {
                     DecompositionFactory_FDRM.lu(origin.numRows, origin.numCols).apply { decompose(origin.copy()) }
                 }
 
-                override val l: Matrix<Float> by lazy {
+                override val l: Matrix<Float32>
-                    lup.getLower(null).wrapMatrix().withFeature(LFeature)
+                    get() = lup.getLower(null).wrapMatrix().withAttribute(LowerTriangular)
-                }
 
-                override val u: Matrix<Float> by lazy {
-                    lup.getUpper(null).wrapMatrix().withFeature(UFeature)
-                }
 
-                override val p: Matrix<Float> by lazy { lup.getRowPivot(null).wrapMatrix() }
+                override val u: Matrix<Float32>
+                    get() = lup.getUpper(null).wrapMatrix().withAttribute(UpperTriangular)
+                override val pivot: IntBuffer get() = lup.getRowPivotV(null).asBuffer()
             }
 
             else -> null
-        }?.let{
-            type.cast(it)
         }
+
+        @Suppress("UNCHECKED_CAST")
+        return raw as V?
     }
 
     /**
@@ -558,7 +529,6 @@ public object EjmlLinearSpaceFDRM : EjmlLinearSpace<Float, Float32Field, FMatrix
 }
 
 
-
 /**
  * [EjmlLinearSpace] implementation based on [CommonOps_DSCC], [DecompositionFactory_DSCC] operations and
  * [DMatrixSparseCSC] matrices.
@@ -705,64 +675,52 @@ public object EjmlLinearSpaceDSCC : EjmlLinearSpace<Double, Float64Field, DMatri
 
     override fun Double.times(v: Point<Double>): EjmlDoubleVector<DMatrixSparseCSC> = v * this
 
-    @UnstableKMathAPI
+    override fun <V, A : StructureAttribute<V>> computeAttribute(structure: Structure2D<Double>, attribute: A): V? {
-    override fun <F : StructureFeature> computeFeature(structure: Matrix<Double>, type: KClass<out F>): F? {
-        structure.getFeature(type)?.let { return it }
         val origin = structure.toEjml().origin
 
-        return when (type) {
+        val raw: Any? = when (attribute) {
-            QRDecompositionFeature::class -> object : QRDecompositionFeature<Double> {
+            Inverted -> {
-                private val qr by lazy {
+                val res = DMatrixRMaj(origin.numRows,origin.numCols)
-                    DecompositionFactory_DSCC.qr(FillReducing.NONE).apply { decompose(origin.copy()) }
+                CommonOps_DSCC.invert(origin,res)
-                }
+                res.wrapMatrix()
-
-                override val q: Matrix<Double> by lazy {
-                    qr.getQ(null, false).wrapMatrix().withFeature(OrthogonalFeature)
-                }
-
-                override val r: Matrix<Double> by lazy { qr.getR(null, false).wrapMatrix().withFeature(UFeature) }
             }
 
-            CholeskyDecompositionFeature::class -> object : CholeskyDecompositionFeature<Double> {
+            Determinant -> CommonOps_DSCC.det(origin)
+
+            QR -> object : QRDecomposition<Double> {
+                val ejmlQr by lazy { DecompositionFactory_DSCC.qr(FillReducing.NONE).apply { decompose(origin.copy()) } }
+                override val q: Matrix<Double> get() = ejmlQr.getQ(null, false).wrapMatrix()
+                override val r: Matrix<Double> get() = ejmlQr.getR(null, false).wrapMatrix()
+            }
+
+            Cholesky -> object : CholeskyDecomposition<Double> {
                 override val l: Matrix<Double> by lazy {
                     val cholesky =
                         DecompositionFactory_DSCC.cholesky().apply { decompose(origin.copy()) }
 
-                    (cholesky.getT(null) as DMatrix).wrapMatrix().withFeature(LFeature)
+                    (cholesky.getT(null) as DMatrix).wrapMatrix().withAttribute(LowerTriangular)
                 }
             }
 
-            LUDecompositionFeature::class, DeterminantFeature::class, InverseMatrixFeature::class -> object :
+            LUP -> object : LupDecomposition<Double> {
-                LUDecompositionFeature<Double>, DeterminantFeature<Double>, InverseMatrixFeature<Double> {
+                private val lup by lazy {
-                private val lu by lazy {
                     DecompositionFactory_DSCC.lu(FillReducing.NONE).apply { decompose(origin.copy()) }
                 }
 
-                override val l: Matrix<Double> by lazy {
+                override val l: Matrix<Double>
-                    lu.getLower(null).wrapMatrix().withFeature(LFeature)
+                    get() = lup.getLower(null).wrapMatrix().withAttribute(LowerTriangular)
-                }
 
-                override val u: Matrix<Double> by lazy {
-                    lu.getUpper(null).wrapMatrix().withFeature(UFeature)
-                }
 
-                override val inverse: Matrix<Double> by lazy {
+                override val u: Matrix<Double>
-                    var a = origin
+                    get() = lup.getUpper(null).wrapMatrix().withAttribute(UpperTriangular)
-                    val inverse = DMatrixRMaj(1, 1)
+                override val pivot: IntBuffer get() = lup.getRowPivotV(null).asBuffer()
-                    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)
         }
+
+        @Suppress("UNCHECKED_CAST")
+        return raw as V?
     }
 
     /**
@@ -793,7 +751,6 @@ public object EjmlLinearSpaceDSCC : EjmlLinearSpace<Double, Float64Field, DMatri
 }
 
 
-
 /**
  * [EjmlLinearSpace] implementation based on [CommonOps_FSCC], [DecompositionFactory_FSCC] operations and
  * [FMatrixSparseCSC] matrices.
@@ -939,65 +896,52 @@ public object EjmlLinearSpaceFSCC : EjmlLinearSpace<Float, Float32Field, FMatrix
     }
 
     override fun Float.times(v: Point<Float>): EjmlFloatVector<FMatrixSparseCSC> = v * this
-
+    override fun <V, A : StructureAttribute<V>> computeAttribute(structure: Structure2D<Float32>, attribute: A): V? {
-    @UnstableKMathAPI
-    override fun <F : StructureFeature> computeFeature(structure: Matrix<Float>, type: KClass<out F>): F? {
-        structure.getFeature(type)?.let { return it }
         val origin = structure.toEjml().origin
 
-        return when (type) {
+        val raw: Any? = when (attribute) {
-            QRDecompositionFeature::class -> object : QRDecompositionFeature<Float> {
+            Inverted -> {
-                private val qr by lazy {
+                val res = FMatrixRMaj(origin.numRows,origin.numCols)
-                    DecompositionFactory_FSCC.qr(FillReducing.NONE).apply { decompose(origin.copy()) }
+                CommonOps_FSCC.invert(origin,res)
-                }
+                res.wrapMatrix()
-
-                override val q: Matrix<Float> by lazy {
-                    qr.getQ(null, false).wrapMatrix().withFeature(OrthogonalFeature)
-                }
-
-                override val r: Matrix<Float> by lazy { qr.getR(null, false).wrapMatrix().withFeature(UFeature) }
             }
 
-            CholeskyDecompositionFeature::class -> object : CholeskyDecompositionFeature<Float> {
+            Determinant -> CommonOps_FSCC.det(origin)
-                override val l: Matrix<Float> by lazy {
+
+            QR -> object : QRDecomposition<Float32> {
+                val ejmlQr by lazy { DecompositionFactory_FSCC.qr(FillReducing.NONE).apply { decompose(origin.copy()) } }
+                override val q: Matrix<Float32> get() = ejmlQr.getQ(null, false).wrapMatrix()
+                override val r: Matrix<Float32> get() = ejmlQr.getR(null, false).wrapMatrix()
+            }
+
+            Cholesky -> object : CholeskyDecomposition<Float32> {
+                override val l: Matrix<Float32> by lazy {
                     val cholesky =
                         DecompositionFactory_FSCC.cholesky().apply { decompose(origin.copy()) }
 
-                    (cholesky.getT(null) as FMatrix).wrapMatrix().withFeature(LFeature)
+                    (cholesky.getT(null) as FMatrix).wrapMatrix().withAttribute(LowerTriangular)
                 }
             }
 
-            LUDecompositionFeature::class, DeterminantFeature::class, InverseMatrixFeature::class -> object :
+            LUP -> object : LupDecomposition<Float32> {
-                LUDecompositionFeature<Float>, DeterminantFeature<Float>, InverseMatrixFeature<Float> {
+                private val lup by lazy {
-                private val lu by lazy {
                     DecompositionFactory_FSCC.lu(FillReducing.NONE).apply { decompose(origin.copy()) }
                 }
 
-                override val l: Matrix<Float> by lazy {
+                override val l: Matrix<Float32>
-                    lu.getLower(null).wrapMatrix().withFeature(LFeature)
+                    get() = lup.getLower(null).wrapMatrix().withAttribute(LowerTriangular)
-                }
 
-                override val u: Matrix<Float> by lazy {
-                    lu.getUpper(null).wrapMatrix().withFeature(UFeature)
-                }
 
-                override val inverse: Matrix<Float> by lazy {
+                override val u: Matrix<Float32>
-                    var a = origin
+                    get() = lup.getUpper(null).wrapMatrix().withAttribute(UpperTriangular)
-                    val inverse = FMatrixRMaj(1, 1)
+                override val pivot: IntBuffer get() = lup.getRowPivotV(null).asBuffer()
-                    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)
         }
+
+        @Suppress("UNCHECKED_CAST")
+        return raw as V?
     }
 
     /**

kmath-ejml/src/test/kotlin/space/kscience/kmath/ejml/EjmlMatrixTest.kt

@@ -58,19 +58,19 @@ internal class EjmlMatrixTest {
 
     @OptIn(UnstableKMathAPI::class)
     @Test
-    fun features() {
+    fun features() = EjmlLinearSpaceDDRM {
         val m = randomMatrix
         val w = EjmlDoubleMatrix(m)
-        val det: Determinant<Double> = EjmlLinearSpaceDDRM.attributeForOrNull(w) ?: fail()
+        val det: Double = w.getOrComputeAttribute(Determinant) ?: fail()
-        assertEquals(CommonOps_DDRM.det(m), det.determinant)
+        assertEquals(CommonOps_DDRM.det(m), det)
-        val lup: LupDecompositionAttribute<Double> = EjmlLinearSpaceDDRM.attributeForOrNull(w) ?: fail()
+        val lup: LupDecomposition<Double> = w.getOrComputeAttribute(LUP) ?: fail()
 
         val ludecompositionF64 = DecompositionFactory_DDRM.lu(m.numRows, m.numCols)
             .also { it.decompose(m.copy()) }
 
         assertMatrixEquals(EjmlDoubleMatrix(ludecompositionF64.getLower(null)), lup.l)
         assertMatrixEquals(EjmlDoubleMatrix(ludecompositionF64.getUpper(null)), lup.u)
-        assertMatrixEquals(EjmlDoubleMatrix(ludecompositionF64.getRowPivot(null)), lup.p)
+        assertMatrixEquals(EjmlDoubleMatrix(ludecompositionF64.getRowPivot(null)), lup.pivotMatrix(this))
     }
 
     @Test

kmath-for-real/build.gradle.kts

@@ -6,6 +6,7 @@ kscience {
     jvm()
     js()
     native()
+    wasm()
 
     dependencies {
         api(projects.kmathCore)

kmath-for-real/src/commonTest/kotlin/space/kscience/kmath/real/DoubleVectorTest.kt

@@ -7,7 +7,7 @@ package space.kscience.kmath.real
 
 import space.kscience.kmath.linear.asMatrix
 import space.kscience.kmath.linear.linearSpace
-import space.kscience.kmath.linear.transpose
+import space.kscience.kmath.linear.transposed
 import space.kscience.kmath.operations.algebra
 import space.kscience.kmath.structures.Float64Buffer
 import kotlin.test.Test
@@ -34,7 +34,7 @@ internal class DoubleVectorTest {
         val vector1 = Float64Buffer(5) { it.toDouble() }
         val vector2 = Float64Buffer(5) { 5 - it.toDouble() }
         val matrix1 = vector1.asMatrix()
-        val matrix2 = vector2.asMatrix().transpose()
+        val matrix2 = vector2.asMatrix().transposed()
         val product = matrix1 dot matrix2
         assertEquals(5.0, product[1, 0])
         assertEquals(6.0, product[2, 2])

kmath-functions/build.gradle.kts

@@ -6,7 +6,6 @@ kscience{
     jvm()
     js()
     native()
-
     wasm()
 
     dependencies {

kmath-functions/src/commonMain/kotlin/space/kscience/kmath/integration/GaussIntegrator.kt

@@ -92,7 +92,7 @@ public inline fun <reified T : Any> GaussIntegrator<T>.integrate(
     range: ClosedRange<Double>,
     order: Int = 10,
     intervals: Int = 10,
-    attributesBuilder: AttributesBuilder<UnivariateIntegrand<T>>.() -> Unit,
+    attributesBuilder: AttributesBuilder<UnivariateIntegrand<T>>.() -> Unit = {},
     noinline function: (Double) -> T,
 ): UnivariateIntegrand<T> {
     require(range.endInclusive > range.start) { "The range upper bound should be higher than lower bound" }

kmath-functions/src/commonMain/kotlin/space/kscience/kmath/interpolation/SplineInterpolator.kt

@@ -23,7 +23,7 @@ import space.kscience.kmath.structures.MutableBufferFactory
  */
 public class SplineInterpolator<T : Comparable<T>>(
     override val algebra: Field<T>,
-    public val bufferFactory: MutableBufferFactory<T>,
+    public val bufferFactory: MutableBufferFactory<T> = algebra.bufferFactory,
 ) : PolynomialInterpolator<T> {
     //TODO possibly optimize zeroed buffers
 

kmath-functions/src/commonTest/kotlin/space/kscience/kmath/functions/testUtils/IntModulo.kt

@@ -9,6 +9,7 @@ package space.kscience.kmath.functions.testUtils
 
 import space.kscience.kmath.operations.Ring
 import space.kscience.kmath.operations.ScaleOperations
+import space.kscience.kmath.structures.MutableBufferFactory
 
 
 class IntModulo {
@@ -109,15 +110,17 @@ class IntModulo {
 }
 
 @Suppress("EXTENSION_SHADOWED_BY_MEMBER", "OVERRIDE_BY_INLINE")
-class IntModuloRing : Ring<IntModulo>, ScaleOperations<IntModulo> {
+class IntModuloRing(modulus: Int) : Ring<IntModulo>, ScaleOperations<IntModulo> {
 
     val modulus: Int
 
-    constructor(modulus: Int) {
+    init {
         require(modulus != 0) { "modulus can not be zero" }
         this.modulus = if (modulus < 0) -modulus else modulus
     }
 
+    override val bufferFactory: MutableBufferFactory<IntModulo> = MutableBufferFactory()
+
     override inline val zero: IntModulo get() = IntModulo(0, modulus, toCheckInput = false)
     override inline val one: IntModulo get() = IntModulo(1, modulus, toCheckInput = false)
 

kmath-functions/src/commonTest/kotlin/space/kscience/kmath/functions/testUtils/Rational.kt

@@ -10,6 +10,7 @@ package space.kscience.kmath.functions.testUtils
 import space.kscience.kmath.UnstableKMathAPI
 import space.kscience.kmath.operations.Field
 import space.kscience.kmath.operations.NumbersAddOps
+import space.kscience.kmath.structures.MutableBufferFactory
 
 @Suppress("NAME_SHADOWING")
 class Rational {
@@ -159,6 +160,7 @@ class Rational {
 @Suppress("EXTENSION_SHADOWED_BY_MEMBER", "OVERRIDE_BY_INLINE")
 @OptIn(UnstableKMathAPI::class)
 object RationalField : Field<Rational>, NumbersAddOps<Rational> {
+    override val bufferFactory: MutableBufferFactory<Rational> = MutableBufferFactory()
     override inline val zero: Rational get() = Rational.ZERO
     override inline val one: Rational get() = Rational.ONE
 

kmath-functions/src/commonTest/kotlin/space/kscience/kmath/integration/SplineIntegralTest.kt

@@ -18,15 +18,15 @@ import kotlin.test.assertEquals
 class SplineIntegralTest {
 
     @Test
-    fun integratePolynomial(){
+    fun integratePolynomial() {
         val polynomial = Polynomial(1.0, 2.0, 3.0)
-        val integral = polynomial.integrate(Float64Field,1.0..2.0)
+        val integral = polynomial.integrate(Float64Field, 1.0..2.0)
         assertEquals(11.0, integral, 0.001)
     }
 
     @Test
     fun gaussSin() {
-        val res = Float64Field.splineIntegrator.integrate(0.0..2 * PI, IntegrandMaxCalls(5)) { x ->
+        val res = Float64Field.splineIntegrator.integrate(0.0..2 * PI, { IntegrandMaxCalls(5) }) { x ->
             sin(x)
         }
         assertEquals(0.0, res.value, 1e-2)
@@ -34,8 +34,8 @@ class SplineIntegralTest {
 
     @Test
     fun gaussUniform() {
-        val res = Float64Field.splineIntegrator.integrate(35.0..100.0, IntegrandMaxCalls(20)) { x ->
+        val res = Float64Field.splineIntegrator.integrate(35.0..100.0, { IntegrandMaxCalls(20) }) { x ->
-            if(x in 30.0..50.0){
+            if (x in 30.0..50.0) {
                 1.0
             } else {
                 0.0

kmath-geometry/src/commonMain/kotlin/space/kscience/kmath/geometry/Line.kt

@@ -7,8 +7,8 @@ package space.kscience.kmath.geometry
 
 import kotlinx.serialization.SerialName
 import kotlinx.serialization.Serializable
-import space.kscience.kmath.geometry.euclidean2d.DoubleVector2D
+import space.kscience.kmath.geometry.euclidean2d.Float64Vector2D
-import space.kscience.kmath.geometry.euclidean3d.DoubleVector3D
+import space.kscience.kmath.geometry.euclidean3d.Float64Vector3D
 
 /**
  * A line formed by [start] vector of start and a [direction] vector. Direction vector is not necessarily normalized,
@@ -25,8 +25,8 @@ private data class LineImpl<out V : Any>(override val start: V, override val dir
 
 public fun <V : Any> Line(base: V, direction: V): Line<V> = LineImpl(base, direction)
 
-public typealias Line2D = Line<DoubleVector2D>
+public typealias Line2D = Line<Float64Vector2D>
-public typealias Line3D = Line<DoubleVector3D>
+public typealias Line3D = Line<Float64Vector3D>
 
 /**
  * A directed line segment between [begin] and [end]
@@ -49,5 +49,5 @@ public fun <V : Any> LineSegment<V>.line(algebra: GeometrySpace<V, *>): Line<V>
     Line(begin, end - begin)
 }
 
-public typealias LineSegment2D = LineSegment<DoubleVector2D>
+public typealias LineSegment2D = LineSegment<Float64Vector2D>
-public typealias LineSegment3D = LineSegment<DoubleVector3D>
+public typealias LineSegment3D = LineSegment<Float64Vector3D>

kmath-geometry/src/commonMain/kotlin/space/kscience/kmath/geometry/Vector3D.kt

@@ -5,6 +5,7 @@
 
 package space.kscience.kmath.geometry
 
+import space.kscience.attributes.SafeType
 import space.kscience.kmath.linear.Point
 import space.kscience.kmath.structures.Buffer
 
@@ -33,6 +34,9 @@ public fun <T> Buffer<T>.asVector3D(): Vector3D<T> = object : Vector3D<T> {
         require(this@asVector3D.size == 3) { "Buffer of size 3 is required for Vector3D" }
     }
 
+    override val type: SafeType<T> = this@asVector3D.type
+
+
     override val x: T get() = this@asVector3D[0]
     override val y: T get() = this@asVector3D[1]
     override val z: T get() = this@asVector3D[2]

kmath-geometry/src/commonMain/kotlin/space/kscience/kmath/geometry/angles.kt

@@ -11,6 +11,8 @@ import kotlinx.serialization.builtins.serializer
 import kotlinx.serialization.descriptors.SerialDescriptor
 import kotlinx.serialization.encoding.Decoder
 import kotlinx.serialization.encoding.Encoder
+import space.kscience.kmath.operations.Group
+import space.kscience.kmath.structures.MutableBufferFactory
 import kotlin.jvm.JvmInline
 import kotlin.math.PI
 import kotlin.math.floor
@@ -28,11 +30,19 @@ public sealed interface Angle : Comparable<Angle> {
     public operator fun div(other: Angle): Double
     public operator fun unaryMinus(): Angle
 
-    public companion object {
+    public companion object: Group<Angle> {
-        public val zero: Radians = Radians(0.0)
+        override val zero: Radians = Radians(0.0)
         public val pi: Radians = Radians(PI)
         public val piTimes2: Radians = Radians(PI * 2)
         public val piDiv2: Radians = Radians(PI / 2)
+
+
+        override fun add(left: Angle, right: Angle): Angle = left + right
+
+        @Suppress("EXTENSION_SHADOWED_BY_MEMBER")
+        override fun Angle.unaryMinus(): Angle = -this
+
+        override val bufferFactory: MutableBufferFactory<Angle> = MutableBufferFactory()
     }
 }
 
@@ -43,7 +53,7 @@ public object AngleSerializer : KSerializer<Angle> {
     override fun deserialize(decoder: Decoder): Angle = decoder.decodeDouble().degrees
 
     override fun serialize(encoder: Encoder, value: Angle) {
-        encoder.encodeDouble(value.degrees)
+        encoder.encodeDouble(value.toDegrees().value)
     }
 }
 
@@ -56,16 +66,16 @@ public value class Radians(public val value: Double) : Angle {
     override fun toRadians(): Radians = this
     override fun toDegrees(): Degrees = Degrees(value * 180 / PI)
 
-    public override fun plus(other: Angle): Radians = Radians(value + other.radians)
+    public override fun plus(other: Angle): Radians = Radians(value + other.toRadians().value)
-    public override fun minus(other: Angle): Radians = Radians(value - other.radians)
+    public override fun minus(other: Angle): Radians = Radians(value - other.toRadians().value)
 
     public override fun times(other: Number): Radians = Radians(value * other.toDouble())
     public override fun div(other: Number): Radians = Radians(value / other.toDouble())
-    override fun div(other: Angle): Double = value / other.radians
+    override fun div(other: Angle): Double = value / other.toRadians().value
 
     public override fun unaryMinus(): Radians = Radians(-value)
 
-    override fun compareTo(other: Angle): Int = value.compareTo(other.radians)
+    override fun compareTo(other: Angle): Int = value.compareTo(other.toRadians().value)
 }
 
 public fun sin(angle: Angle): Double = kotlin.math.sin(angle.toRadians().value)
@@ -85,16 +95,16 @@ public value class Degrees(public val value: Double) : Angle {
     override fun toRadians(): Radians = Radians(value * PI / 180)
     override fun toDegrees(): Degrees = this
 
-    public override fun plus(other: Angle): Degrees = Degrees(value + other.degrees)
+    public override fun plus(other: Angle): Degrees = Degrees(value + other.toDegrees().value)
-    public override fun minus(other: Angle): Degrees = Degrees(value - other.degrees)
+    public override fun minus(other: Angle): Degrees = Degrees(value - other.toDegrees().value)
 
     public override fun times(other: Number): Degrees = Degrees(value * other.toDouble())
     public override fun div(other: Number): Degrees = Degrees(value / other.toDouble())
-    override fun div(other: Angle): Double = value / other.degrees
+    override fun div(other: Angle): Double = value / other.toDegrees().value
 
     public override fun unaryMinus(): Degrees = Degrees(-value)
 
-    override fun compareTo(other: Angle): Int = value.compareTo(other.degrees)
+    override fun compareTo(other: Angle): Int = value.compareTo(other.toDegrees().value)
 }
 
 public val Number.degrees: Degrees get() = Degrees(toDouble())
@@ -106,6 +116,6 @@ public val Angle.degrees: Double get() = toDegrees().value
  * Normalized angle 2 PI range symmetric around [center]. By default, uses (0, 2PI) range.
  */
 public fun Angle.normalized(center: Angle = Angle.pi): Angle =
-    this - Angle.piTimes2 * floor((radians + PI - center.radians) / PI / 2)
+    this - Angle.piTimes2 * floor((toRadians().value + PI - center.toRadians().value) / PI / 2)
 
 public fun abs(angle: Angle): Angle = if (angle < Angle.zero) -angle else angle

kmath-geometry/src/commonMain/kotlin/space/kscience/kmath/geometry/euclidean2d/Float32Space2D.kt

@@ -11,19 +11,22 @@ import kotlinx.serialization.Serializable
 import kotlinx.serialization.descriptors.SerialDescriptor
 import kotlinx.serialization.encoding.Decoder
 import kotlinx.serialization.encoding.Encoder
+import space.kscience.attributes.SafeType
 import space.kscience.kmath.geometry.GeometrySpace
 import space.kscience.kmath.geometry.Vector2D
 import space.kscience.kmath.operations.Float32Field
 import space.kscience.kmath.structures.Float32
+import space.kscience.kmath.structures.MutableBufferFactory
 import kotlin.math.pow
 import kotlin.math.sqrt
 
 @Serializable(Float32Space2D.VectorSerializer::class)
-public interface Float32Vector2D : Vector2D<Float>
+public interface Float32Vector2D : Vector2D<Float32>{
+    override val type: SafeType<Float32> get() = Float32Field.type
+}
 
 
 public object Float32Space2D : GeometrySpace<Float32Vector2D, Float32> {
-
     @Serializable
     @SerialName("Float32Vector2D")
     private data class Vector2DImpl(
@@ -72,6 +75,8 @@ public object Float32Space2D : GeometrySpace<Float32Vector2D, Float32> {
     public val yAxis: Float32Vector2D = vector(0.0, 1.0)
 
     override val defaultPrecision: Float32 = 1e-3f
+
+    override val bufferFactory: MutableBufferFactory<Float32Vector2D> = MutableBufferFactory()
 }
 
 public fun Float32Vector2D(x: Number, y: Number): Float32Vector2D = Float32Space2D.vector(x, y)

kmath-geometry/src/commonMain/kotlin/space/kscience/kmath/geometry/euclidean2d/Float64Space2D.kt

@@ -11,66 +11,77 @@ import kotlinx.serialization.Serializable
 import kotlinx.serialization.descriptors.SerialDescriptor
 import kotlinx.serialization.encoding.Decoder
 import kotlinx.serialization.encoding.Encoder
+import space.kscience.attributes.SafeType
 import space.kscience.kmath.geometry.GeometrySpace
 import space.kscience.kmath.geometry.Vector2D
-import space.kscience.kmath.linear.Float64LinearSpace
 import space.kscience.kmath.operations.Float64Field
 import space.kscience.kmath.operations.ScaleOperations
+import space.kscience.kmath.structures.Float64
+import space.kscience.kmath.structures.MutableBufferFactory
 import kotlin.math.pow
 import kotlin.math.sqrt
 
 
-public typealias DoubleVector2D = Vector2D<Double>
+@Serializable(Float64Space2D.VectorSerializer::class)
-public typealias Float64Vector2D = Vector2D<Double>
+public interface Float64Vector2D : Vector2D<Float64> {
+    override val type: SafeType<Float64> get() = Float64Field.type
+}
+
+@Deprecated("Use Float64Vector2D", ReplaceWith("Float64Vector2D"))
+public typealias DoubleVector2D = Float64Vector2D
 
-public val Vector2D<Double>.r: Double get() = Float64Space2D.norm(this)
 
 
 /**
  * 2D Euclidean space
  */
-public object Float64Space2D : GeometrySpace<DoubleVector2D>, ScaleOperations<DoubleVector2D> {
+public object Float64Space2D : GeometrySpace<Float64Vector2D, Float64>, ScaleOperations<Float64Vector2D> {
+
 
     @Serializable
     @SerialName("Float64Vector2D")
     private data class Vector2DImpl(
         override val x: Double,
         override val y: Double,
-    ) : DoubleVector2D
+    ) : Float64Vector2D
 
-    public object VectorSerializer : KSerializer<DoubleVector2D> {
+    public object VectorSerializer : KSerializer<Float64Vector2D> {
         private val proxySerializer = Vector2DImpl.serializer()
         override val descriptor: SerialDescriptor get() = proxySerializer.descriptor
 
-        override fun deserialize(decoder: Decoder): DoubleVector2D = decoder.decodeSerializableValue(proxySerializer)
+        override fun deserialize(decoder: Decoder): Float64Vector2D = decoder.decodeSerializableValue(proxySerializer)
 
-        override fun serialize(encoder: Encoder, value: DoubleVector2D) {
+        override fun serialize(encoder: Encoder, value: Float64Vector2D) {
             val vector = value as? Vector2DImpl ?: Vector2DImpl(value.x, value.y)
             encoder.encodeSerializableValue(proxySerializer, vector)
         }
     }
 
-    public fun vector(x: Number, y: Number): DoubleVector2D = Vector2DImpl(x.toDouble(), y.toDouble())
+    public fun vector(x: Number, y: Number): Float64Vector2D = Vector2DImpl(x.toDouble(), y.toDouble())
 
-    override val zero: DoubleVector2D by lazy { vector(0.0, 0.0) }
+    override val zero: Float64Vector2D by lazy { vector(0.0, 0.0) }
 
-    override fun norm(arg: DoubleVector2D): Double = sqrt(arg.x.pow(2) + arg.y.pow(2))
+    override fun norm(arg: Float64Vector2D): Double = sqrt(arg.x.pow(2) + arg.y.pow(2))
 
-    override fun DoubleVector2D.unaryMinus(): DoubleVector2D = vector(-x, -y)
+    override fun Float64Vector2D.unaryMinus(): Float64Vector2D = vector(-x, -y)
 
-    override fun DoubleVector2D.distanceTo(other: DoubleVector2D): Double = norm(this - other)
+    override fun Float64Vector2D.distanceTo(other: Float64Vector2D): Double = norm(this - other)
-    override fun add(left: DoubleVector2D, right: DoubleVector2D): DoubleVector2D =
+    override fun add(left: Float64Vector2D, right: Float64Vector2D): Float64Vector2D =
         vector(left.x + right.x, left.y + right.y)
 
-    override fun scale(a: DoubleVector2D, value: Double): DoubleVector2D = vector(a.x * value, a.y * value)
+    override fun scale(a: Float64Vector2D, value: Double): Float64Vector2D = vector(a.x * value, a.y * value)
-    override fun DoubleVector2D.dot(other: DoubleVector2D): Double = x * other.x + y * other.y
+    override fun Float64Vector2D.dot(other: Float64Vector2D): Double = x * other.x + y * other.y
 
-    public val xAxis: DoubleVector2D = vector(1.0, 0.0)
+    public val xAxis: Float64Vector2D = vector(1.0, 0.0)
-    public val yAxis: DoubleVector2D = vector(0.0, 1.0)
+    public val yAxis: Float64Vector2D = vector(0.0, 1.0)
 
     override val defaultPrecision: Double = 1e-6
+
+    override val bufferFactory: MutableBufferFactory<Float64Vector2D> = MutableBufferFactory()
 }
 
 public fun Float64Vector2D(x: Number, y: Number): Float64Vector2D = Float64Space2D.vector(x, y)
 
+public val Float64Vector2D.r: Float64 get() = Float64Space2D.norm(this)
+
 public val Float64Field.euclidean2D: Float64Space2D get() = Float64Space2D

kmath-geometry/src/commonMain/kotlin/space/kscience/kmath/geometry/euclidean3d/Float32Space3D.kt

@@ -11,15 +11,19 @@ import kotlinx.serialization.Serializable
 import kotlinx.serialization.descriptors.SerialDescriptor
 import kotlinx.serialization.encoding.Decoder
 import kotlinx.serialization.encoding.Encoder
+import space.kscience.attributes.SafeType
 import space.kscience.kmath.geometry.GeometrySpace
 import space.kscience.kmath.geometry.Vector3D
 import space.kscience.kmath.operations.Float32Field
 import space.kscience.kmath.structures.Float32
+import space.kscience.kmath.structures.MutableBufferFactory
 import kotlin.math.pow
 import kotlin.math.sqrt
 
 @Serializable(Float32Space3D.VectorSerializer::class)
-public interface Float32Vector3D : Vector3D<Float>
+public interface Float32Vector3D : Vector3D<Float>{
+    override val type: SafeType<Float32> get() = Float32Field.type
+}
 
 
 public object Float32Space3D : GeometrySpace<Float32Vector3D, Float32> {
@@ -101,6 +105,8 @@ public object Float32Space3D : GeometrySpace<Float32Vector3D, Float32> {
     public val zAxis: Float32Vector3D = vector(0.0, 0.0, 1.0)
 
     override val defaultPrecision: Float32 = 1e-3f
+
+    override val bufferFactory: MutableBufferFactory<Float32Vector3D> = MutableBufferFactory()
 }
 
 public fun Float32Vector3D(x: Number, y: Number, z: Number): Float32Vector3D = Float32Space3D.vector(x, y, z)

kmath-geometry/src/commonMain/kotlin/space/kscience/kmath/geometry/euclidean3d/Float64Space3D.kt

@@ -11,10 +11,13 @@ import kotlinx.serialization.Serializable
 import kotlinx.serialization.descriptors.SerialDescriptor
 import kotlinx.serialization.encoding.Decoder
 import kotlinx.serialization.encoding.Encoder
+import space.kscience.attributes.SafeType
 import space.kscience.kmath.geometry.GeometrySpace
 import space.kscience.kmath.geometry.Vector3D
 import space.kscience.kmath.linear.Float64LinearSpace
 import space.kscience.kmath.operations.Float64Field
+import space.kscience.kmath.structures.Float64
+import space.kscience.kmath.structures.MutableBufferFactory
 import kotlin.math.pow
 import kotlin.math.sqrt
 
@@ -31,12 +34,16 @@ internal fun leviCivita(i: Int, j: Int, k: Int): Int = when {
     else -> 0
 }
 
-public typealias DoubleVector3D = Vector3D<Double>
+@Serializable(Float64Space3D.VectorSerializer::class)
-public typealias Float64Vector3D = Vector3D<Double>
+public interface Float64Vector3D : Vector3D<Float64> {
+    override val type: SafeType<Float64> get() = Float64Field.type
+}
 
-public val DoubleVector3D.r: Double get() = Float64Space3D.norm(this)
+@Deprecated("Use Float64Vector3D", ReplaceWith("Float64Vector3D"))
+public typealias DoubleVector3D = Float64Vector3D
 
-public object Float64Space3D : GeometrySpace<DoubleVector3D, Double>{
+
+public object Float64Space3D : GeometrySpace<Vector3D<Float64>, Double> {
 
     public val linearSpace: Float64LinearSpace = Float64LinearSpace
 
@@ -46,52 +53,52 @@ public object Float64Space3D : GeometrySpace<DoubleVector3D, Double>{
         override val x: Double,
         override val y: Double,
         override val z: Double,
-    ) : DoubleVector3D
+    ) : Float64Vector3D
 
-    public object VectorSerializer : KSerializer<DoubleVector3D> {
+    public object VectorSerializer : KSerializer<Float64Vector3D> {
         private val proxySerializer = Vector3DImpl.serializer()
         override val descriptor: SerialDescriptor get() = proxySerializer.descriptor
 
-        override fun deserialize(decoder: Decoder): DoubleVector3D = decoder.decodeSerializableValue(proxySerializer)
+        override fun deserialize(decoder: Decoder): Float64Vector3D = decoder.decodeSerializableValue(proxySerializer)
 
-        override fun serialize(encoder: Encoder, value: DoubleVector3D) {
+        override fun serialize(encoder: Encoder, value: Float64Vector3D) {
             val vector = value as? Vector3DImpl ?: Vector3DImpl(value.x, value.y, value.z)
             encoder.encodeSerializableValue(proxySerializer, vector)
         }
     }
 
-    public fun vector(x: Double, y: Double, z: Double): DoubleVector3D =
+    public fun vector(x: Double, y: Double, z: Double): Float64Vector3D =
         Vector3DImpl(x, y, z)
 
-    public fun vector(x: Number, y: Number, z: Number): DoubleVector3D =
+    public fun vector(x: Number, y: Number, z: Number): Float64Vector3D =
         vector(x.toDouble(), y.toDouble(), z.toDouble())
 
-    override val zero: DoubleVector3D by lazy { vector(0.0, 0.0, 0.0) }
+    override val zero: Float64Vector3D by lazy { vector(0.0, 0.0, 0.0) }
 
-    override fun norm(arg: DoubleVector3D): Double = sqrt(arg.x.pow(2) + arg.y.pow(2) + arg.z.pow(2))
+    override fun norm(arg: Vector3D<Float64>): Double = sqrt(arg.x.pow(2) + arg.y.pow(2) + arg.z.pow(2))
 
-    public fun DoubleVector3D.norm(): Double = norm(this)
+    public fun Vector3D<Float64>.norm(): Double = norm(this)
 
-    override fun DoubleVector3D.unaryMinus(): DoubleVector3D = vector(-x, -y, -z)
+    override fun Vector3D<Float64>.unaryMinus(): Float64Vector3D = vector(-x, -y, -z)
 
-    override fun DoubleVector3D.distanceTo(other: DoubleVector3D): Double = (this - other).norm()
+    override fun Vector3D<Float64>.distanceTo(other: Vector3D<Float64>): Double = (this - other).norm()
 
-    override fun add(left: DoubleVector3D, right: DoubleVector3D): DoubleVector3D =
+    override fun add(left: Vector3D<Float64>, right: Vector3D<Float64>): Float64Vector3D =
         vector(left.x + right.x, left.y + right.y, left.z + right.z)
 
-    override fun scale(a: DoubleVector3D, value: Double): DoubleVector3D =
+    override fun scale(a: Vector3D<Float64>, value: Double): Float64Vector3D =
         vector(a.x * value, a.y * value, a.z * value)
 
-    override fun DoubleVector3D.dot(other: DoubleVector3D): Double =
+    override fun Vector3D<Float64>.dot(other: Vector3D<Float64>): Double =
         x * other.x + y * other.y + z * other.z
 
     /**
      * Compute vector product of [first] and [second]. The basis is assumed to be right-handed.
      */
     public fun vectorProduct(
-        first: DoubleVector3D,
+        first: Vector3D<Float64>,
-        second: DoubleVector3D,
+        second: Vector3D<Float64>,
-    ): DoubleVector3D {
+    ): Float64Vector3D {
         var x = 0.0
         var y = 0.0
         var z = 0.0
@@ -110,13 +117,18 @@ public object Float64Space3D : GeometrySpace<DoubleVector3D, Double>{
     /**
      * Vector product with the right basis
      */
-    public infix fun DoubleVector3D.cross(other: DoubleVector3D): Vector3D<Double> = vectorProduct(this, other)
+    public infix fun Vector3D<Float64>.cross(other: Vector3D<Float64>): Vector3D<Double> = vectorProduct(this, other)
 
-    public val xAxis: DoubleVector3D = vector(1.0, 0.0, 0.0)
+    public val xAxis: Float64Vector3D = vector(1.0, 0.0, 0.0)
-    public val yAxis: DoubleVector3D = vector(0.0, 1.0, 0.0)
+    public val yAxis: Float64Vector3D = vector(0.0, 1.0, 0.0)
-    public val zAxis: DoubleVector3D = vector(0.0, 0.0, 1.0)
+    public val zAxis: Float64Vector3D = vector(0.0, 0.0, 1.0)
 
     override val defaultPrecision: Double = 1e-6
+
+    override val bufferFactory: MutableBufferFactory<Vector3D<Float64>> = MutableBufferFactory()
 }
 
 public val Float64Field.euclidean3D: Float64Space3D get() = Float64Space3D
+
+
+public val Float64Vector3D.r: Double get() = Float64Space3D.norm(this)

kmath-geometry/src/commonMain/kotlin/space/kscience/kmath/geometry/euclidean3d/rotations3D.kt

@@ -5,6 +5,7 @@
 
 package space.kscience.kmath.geometry.euclidean3d
 
+import space.kscience.attributes.SafeType
 import space.kscience.kmath.UnstableKMathAPI
 import space.kscience.kmath.complex.*
 import space.kscience.kmath.geometry.*
@@ -13,6 +14,7 @@ import space.kscience.kmath.linear.Matrix
 import space.kscience.kmath.linear.linearSpace
 import space.kscience.kmath.linear.matrix
 import space.kscience.kmath.operations.Float64Field
+import space.kscience.kmath.structures.Float64
 import kotlin.math.*
 
 public operator fun Quaternion.times(other: Quaternion): Quaternion = QuaternionAlgebra.multiply(this, other)
@@ -35,7 +37,7 @@ public infix fun Quaternion.dot(other: Quaternion): Double = w * other.w + x * o
 /**
  * Represent a vector as quaternion with zero a rotation angle.
  */
-internal fun DoubleVector3D.asQuaternion(): Quaternion = Quaternion(0.0, x, y, z)
+internal fun Float64Vector3D.asQuaternion(): Quaternion = Quaternion(0.0, x, y, z)
 
 /**
  * Angle in radians denoted by this quaternion rotation
@@ -45,7 +47,7 @@ public val Quaternion.theta: Radians get() = (kotlin.math.acos(normalized().w) *
 /**
  * Create a normalized Quaternion from rotation angle and rotation vector
  */
-public fun Quaternion.Companion.fromRotation(theta: Angle, vector: DoubleVector3D): Quaternion {
+public fun Quaternion.Companion.fromRotation(theta: Angle, vector: Float64Vector3D): Quaternion {
     val s = sin(theta / 2)
     val c = cos(theta / 2)
     val norm = with(Float64Space3D) { vector.norm() }
@@ -55,9 +57,9 @@ public fun Quaternion.Companion.fromRotation(theta: Angle, vector: DoubleVector3
 /**
  * An axis of quaternion rotation
  */
-public val Quaternion.vector: DoubleVector3D
+public val Quaternion.vector: Float64Vector3D
     get() {
-        return object : DoubleVector3D {
+        return object : Float64Vector3D {
             private val sint2 = sqrt(1 - w * w)
             override val x: Double get() = this@vector.x / sint2
             override val y: Double get() = this@vector.y / sint2
@@ -69,7 +71,7 @@ public val Quaternion.vector: DoubleVector3D
 /**
  * Rotate a vector in a [Float64Space3D] with [quaternion]
  */
-public fun Float64Space3D.rotate(vector: DoubleVector3D, quaternion: Quaternion): DoubleVector3D =
+public fun Float64Space3D.rotate(vector: Float64Vector3D, quaternion: Quaternion): Float64Vector3D =
     with(QuaternionAlgebra) {
         val p = vector.asQuaternion()
         (quaternion * p * quaternion.reciprocal).vector
@@ -80,15 +82,15 @@ public fun Float64Space3D.rotate(vector: DoubleVector3D, quaternion: Quaternion)
  */
 @UnstableKMathAPI
 public fun Float64Space3D.rotate(
-    vector: DoubleVector3D,
+    vector: Float64Vector3D,
     composition: QuaternionAlgebra.() -> Quaternion,
-): DoubleVector3D =
+): Float64Vector3D =
     rotate(vector, QuaternionAlgebra.composition())
 
 /**
  * Rotate a [Float64] vector in 3D space with a rotation matrix
  */
-public fun Float64Space3D.rotate(vector: DoubleVector3D, matrix: Matrix<Double>): DoubleVector3D {
+public fun Float64Space3D.rotate(vector: Float64Vector3D, matrix: Matrix<Double>):  Vector3D<Float64> {
     require(matrix.colNum == 3 && matrix.rowNum == 3) { "Square 3x3 rotation matrix is required" }
     return with(linearSpace) { (matrix dot vector).asVector3D() }
 }
@@ -242,6 +244,8 @@ public fun Quaternion.Companion.fromEuler(
  * A vector consisting of angles
  */
 public data class AngleVector(override val x: Angle, override val y: Angle, override val z: Angle) : Vector3D<Angle> {
+    override val type: SafeType<Angle> get() = Angle.type
+
     public companion object
 }
 

kmath-geometry/src/commonTest/kotlin/space/kscience/kmath/geometry/testUtils.kt

@@ -5,8 +5,7 @@
 
 package space.kscience.kmath.geometry
 
-import space.kscience.kmath.geometry.euclidean2d.DoubleVector2D
+import space.kscience.kmath.structures.Float64
-import space.kscience.kmath.geometry.euclidean3d.DoubleVector3D
 import kotlin.math.abs
 import kotlin.test.assertEquals
 
@@ -27,12 +26,12 @@ fun grid(
     return xs.flatMap { x -> ys.map { y -> x to y } }
 }
 
-fun assertVectorEquals(expected: DoubleVector2D, actual: DoubleVector2D, absoluteTolerance: Double = 1e-3) {
+fun assertVectorEquals(expected:  Vector2D<Float64>, actual:  Vector2D<Float64>, absoluteTolerance: Double = 1e-3) {
     assertEquals(expected.x, actual.x, absoluteTolerance)
     assertEquals(expected.y, actual.y, absoluteTolerance)
 }
 
-fun assertVectorEquals(expected: DoubleVector3D, actual: DoubleVector3D, absoluteTolerance: Double = 1e-6) {
+fun assertVectorEquals(expected: Vector3D<Float64>, actual:  Vector3D<Float64>, absoluteTolerance: Double = 1e-6) {
     assertEquals(expected.x, actual.x, absoluteTolerance)
     assertEquals(expected.y, actual.y, absoluteTolerance)
     assertEquals(expected.z, actual.z, absoluteTolerance)

kmath-histograms/build.gradle.kts

@@ -6,6 +6,8 @@ kscience{
     jvm()
     js()
     native()
+    wasm()
+    useCoroutines()
 }
 
 //apply(plugin = "kotlinx-atomicfu")
@@ -21,7 +23,6 @@ kotlin.sourceSets {
         dependencies {
             implementation(project(":kmath-for-real"))
             implementation(projects.kmath.kmathStat)
-            implementation("org.jetbrains.kotlinx:kotlinx-coroutines-test:1.6.0")
         }
     }
 }

kmath-histograms/src/commonMain/kotlin/space/kscience/kmath/histogram/Histogram.kt

@@ -47,7 +47,7 @@ public interface Histogram<in T : Any, out V, out B : Bin<T, V>> {
     }
 }
 
-public interface HistogramBuilder<in T : Any, V : Any> {
+public interface HistogramBuilder<in T, V> {
 
     /**
      * The default value increment for a bin
@@ -61,9 +61,9 @@ public interface HistogramBuilder<in T : Any, V : Any> {
 
 }
 
-public fun <T : Any> HistogramBuilder<T, *>.put(point: Point<out T>): Unit = putValue(point)
+public fun <T> HistogramBuilder<T, *>.put(point: Point<T>): Unit = putValue(point)
 
-public fun <T : Any> HistogramBuilder<T, *>.put(vararg point: T): Unit = put(point.asBuffer())
+public inline fun <reified T> HistogramBuilder<T, *>.put(vararg point: T): Unit = put(point.asList().asBuffer())
 
 public fun HistogramBuilder<Double, *>.put(vararg point: Number): Unit =
     put(Float64Buffer(point.map { it.toDouble() }.toDoubleArray()))

kmath-histograms/src/commonMain/kotlin/space/kscience/kmath/histogram/UniformHistogram1D.kt

@@ -12,6 +12,7 @@ import space.kscience.kmath.operations.Ring
 import space.kscience.kmath.operations.ScaleOperations
 import space.kscience.kmath.operations.invoke
 import space.kscience.kmath.structures.Buffer
+import space.kscience.kmath.structures.MutableBufferFactory
 import kotlin.math.floor
 
 @OptIn(UnstableKMathAPI::class)
@@ -46,6 +47,8 @@ public class UniformHistogram1DGroup<V : Any, A>(
     public val startPoint: Double = 0.0,
 ) : Group<Histogram1D<Double, V>>, ScaleOperations<Histogram1D<Double, V>> where A : Ring<V>, A : ScaleOperations<V> {
 
+    override val bufferFactory: MutableBufferFactory<Histogram1D<Double, V>> = MutableBufferFactory()
+
     override val zero: UniformHistogram1D<V> = UniformHistogram1D(this, emptyMap())
 
     /**

kmath-histograms/src/commonMain/kotlin/space/kscience/kmath/histogram/UniformHistogramGroupND.kt

@@ -37,6 +37,8 @@ public class UniformHistogramGroupND<V : Any, A : Field<V>>(
         require(!lower.indices.any { upper[it] - lower[it] < 0 }) { "Range for one of axis is not strictly positive" }
     }
 
+    override val bufferFactory: MutableBufferFactory<HistogramND<Double, HyperSquareDomain, V>> = MutableBufferFactory()
+
     public val dimension: Int get() = lower.size
 
     override val shape: ShapeND = ShapeND(IntArray(binNums.size) { binNums[it] + 2 })
@@ -87,7 +89,7 @@ public class UniformHistogramGroupND<V : Any, A : Field<V>>(
         builder: HistogramBuilder<Double, V>.() -> Unit,
     ): HistogramND<Double, HyperSquareDomain, V> {
         val ndCounter: BufferND<ObjectCounter<V>> =
-            StructureND.buffered(shape) { Counter.of(valueAlgebraND.elementAlgebra) }
+            BufferND(shape) { Counter.of(valueAlgebraND.elementAlgebra) }
         val hBuilder = object : HistogramBuilder<Double, V> {
             override val defaultValue: V get() = valueAlgebraND.elementAlgebra.one
 
@@ -97,7 +99,8 @@ public class UniformHistogramGroupND<V : Any, A : Field<V>>(
             }
         }
         hBuilder.apply(builder)
-        val values: BufferND<V> = BufferND(ndCounter.indices, ndCounter.buffer.mapToBuffer(valueBufferFactory) { it.value })
+        val values: BufferND<V> =
+            BufferND(ndCounter.indices, ndCounter.buffer.mapToBuffer(valueBufferFactory) { it.value })
 
         return HistogramND(this, values)
     }
@@ -128,8 +131,7 @@ public fun <V : Any, A : Field<V>> Histogram.Companion.uniformNDFromRanges(
 
 public fun Histogram.Companion.uniformDoubleNDFromRanges(
     vararg ranges: ClosedFloatingPointRange<Double>,
-): UniformHistogramGroupND<Double, Float64Field> =
+): UniformHistogramGroupND<Double, Float64Field> = uniformNDFromRanges(Floa64FieldOpsND, *ranges)
-    uniformNDFromRanges(Floa64FieldOpsND, *ranges, bufferFactory = ::Float64Buffer)
 
 
 /**
@@ -147,21 +149,18 @@ public fun <V : Any, A : Field<V>> Histogram.Companion.uniformNDFromRanges(
     bufferFactory: BufferFactory<V> = valueAlgebraND.elementAlgebra.bufferFactory,
 ): UniformHistogramGroupND<V, A> = UniformHistogramGroupND(
     valueAlgebraND,
-    DoubleBuffer(
+    ranges
-        ranges
+        .map(Pair<ClosedFloatingPointRange<Double>, Int>::first)
-            .map(Pair<ClosedFloatingPointRange<Double>, Int>::first)
+        .map(ClosedFloatingPointRange<Double>::start)
-            .map(ClosedFloatingPointRange<Double>::start)
+        .asBuffer(),
-    ),
+    ranges
-    DoubleBuffer(
+        .map(Pair<ClosedFloatingPointRange<Double>, Int>::first)
-        ranges
+        .map(ClosedFloatingPointRange<Double>::endInclusive)
-            .map(Pair<ClosedFloatingPointRange<Double>, Int>::first)
+        .asBuffer(),
-            .map(ClosedFloatingPointRange<Double>::endInclusive)
-    ),
     ranges.map(Pair<ClosedFloatingPointRange<Double>, Int>::second).toIntArray(),
     valueBufferFactory = bufferFactory
 )
 
 public fun Histogram.Companion.uniformDoubleNDFromRanges(
     vararg ranges: Pair<ClosedFloatingPointRange<Double>, Int>,
-): UniformHistogramGroupND<Double, Float64Field> =
+): UniformHistogramGroupND<Double, Float64Field> = uniformNDFromRanges(Floa64FieldOpsND, *ranges)
-    uniformNDFromRanges(Floa64FieldOpsND, *ranges, bufferFactory = ::Float64Buffer)

kmath-histograms/src/jvmMain/kotlin/space/kscience/kmath/histogram/TreeHistogramGroup.kt

@@ -14,10 +14,7 @@ import space.kscience.kmath.misc.sorted
 import space.kscience.kmath.operations.Group
 import space.kscience.kmath.operations.Ring
 import space.kscience.kmath.operations.ScaleOperations
-import space.kscience.kmath.structures.Buffer
+import space.kscience.kmath.structures.*
-import space.kscience.kmath.structures.first
-import space.kscience.kmath.structures.indices
-import space.kscience.kmath.structures.last
 import java.util.*
 
 private fun <B : ClosedRange<Double>> TreeMap<Double, B>.getBin(value: Double): B? {
@@ -53,6 +50,8 @@ public class TreeHistogramGroup<V : Any, A>(
     @PublishedApi internal val binFactory: (Double) -> DoubleDomain1D,
 ) : Group<TreeHistogram<V>>, ScaleOperations<TreeHistogram<V>> where A : Ring<V>, A : ScaleOperations<V> {
 
+    override val bufferFactory: MutableBufferFactory<TreeHistogram<V>> = MutableBufferFactory()
+
     internal inner class DomainCounter(val domain: DoubleDomain1D, val counter: Counter<V> = Counter.of(valueAlgebra)) :
         ClosedRange<Double> by domain.range
 

kmath-kotlingrad/src/main/kotlin/space/kscience/kmath/kotlingrad/KotlingradExpression.kt

@@ -7,6 +7,7 @@ package space.kscience.kmath.kotlingrad
 
 import ai.hypergraph.kotlingrad.api.SFun
 import ai.hypergraph.kotlingrad.api.SVar
+import space.kscience.attributes.SafeType
 import space.kscience.kmath.expressions.*
 import space.kscience.kmath.operations.NumericAlgebra
 
@@ -25,6 +26,8 @@ public class KotlingradExpression<T : Number, A : NumericAlgebra<T>>(
     public val algebra: A,
     public val mst: MST,
 ) : SpecialDifferentiableExpression<T, KotlingradExpression<T, A>> {
+    override val type: SafeType<T> = algebra.type
+
     override fun invoke(arguments: Map<Symbol, T>): T = mst.interpret(algebra, arguments)
 
     override fun derivativeOrNull(

kmath-nd4j/build.gradle.kts

@@ -12,7 +12,7 @@ dependencies {
 }
 
 readme {
-    maturity = space.kscience.gradle.Maturity.EXPERIMENTAL
+    maturity = space.kscience.gradle.Maturity.DEPRECATED
     propertyByTemplate("artifact", rootProject.file("docs/templates/ARTIFACT-TEMPLATE.md"))
     feature(id = "nd4jarraystructure") { "NDStructure wrapper for INDArray" }
     feature(id = "nd4jarrayrings") { "Rings over Nd4jArrayStructure of Int and Long" }

kmath-nd4j/src/main/kotlin/space/kscience/kmath/nd4j/Nd4jArrayStructure.kt

@@ -6,8 +6,10 @@
 package space.kscience.kmath.nd4j
 
 import org.nd4j.linalg.api.ndarray.INDArray
+import space.kscience.attributes.SafeType
 import space.kscience.kmath.PerformancePitfall
 import space.kscience.kmath.nd.*
+import space.kscience.kmath.operations.Float32Field
 
 /**
  * Represents a [StructureND] wrapping an [INDArray] object.
@@ -31,6 +33,7 @@ public sealed class Nd4jArrayStructure<T> : MutableStructureND<T> {
 }
 
 public data class Nd4jArrayIntStructure(override val ndArray: INDArray) : Nd4jArrayStructure<Int>(), StructureNDOfInt {
+
     override fun elementsIterator(): Iterator<Pair<IntArray, Int>> = ndArray.intIterator()
 
     @OptIn(PerformancePitfall::class)
@@ -47,8 +50,10 @@ public data class Nd4jArrayIntStructure(override val ndArray: INDArray) : Nd4jAr
  */
 public fun INDArray.asIntStructure(): Nd4jArrayIntStructure = Nd4jArrayIntStructure(this)
 
-public data class Nd4jArrayDoubleStructure(override val ndArray: INDArray) : Nd4jArrayStructure<Double>(), StructureNDOfDouble {
+public data class Nd4jArrayDoubleStructure(override val ndArray: INDArray) : Nd4jArrayStructure<Double>(),
+    StructureNDOfDouble {
     override fun elementsIterator(): Iterator<Pair<IntArray, Double>> = ndArray.realIterator()
+
     @OptIn(PerformancePitfall::class)
     override fun get(index: IntArray): Double = ndArray.getDouble(*index)
 
@@ -64,7 +69,11 @@ public data class Nd4jArrayDoubleStructure(override val ndArray: INDArray) : Nd4
 public fun INDArray.asDoubleStructure(): Nd4jArrayStructure<Double> = Nd4jArrayDoubleStructure(this)
 
 public data class Nd4jArrayFloatStructure(override val ndArray: INDArray) : Nd4jArrayStructure<Float>() {
+
+    override val type: SafeType<Float> get() = Float32Field.type
+
     override fun elementsIterator(): Iterator<Pair<IntArray, Float>> = ndArray.floatIterator()
+
     @PerformancePitfall
     override fun get(index: IntArray): Float = ndArray.getFloat(*index)
 

kmath-optimization/build.gradle.kts

@@ -6,6 +6,7 @@ kscience{
     jvm()
     js()
     native()
+    wasm()
 }
 
 kotlin.sourceSets {

kmath-optimization/src/commonMain/kotlin/space/kscience/kmath/optimization/FunctionOptimization.kt

@@ -11,6 +11,10 @@ import space.kscience.kmath.expressions.Symbol
 
 public class OptimizationValue<V>(type: SafeType<V>) : PolymorphicAttribute<V>(type)
 
+public inline fun <reified T> AttributesBuilder<FunctionOptimization<T>>.value(value: T) {
+    set(OptimizationValue(safeTypeOf<T>()), value)
+}
+
 public enum class OptimizationDirection {
     MAXIMIZE,
     MINIMIZE

kmath-optimization/src/commonMain/kotlin/space/kscience/kmath/optimization/OptimizationProblem.kt

@@ -34,10 +34,18 @@ public fun <T> AttributesBuilder<OptimizationProblem<T>>.startAt(startingPoint:
 public class OptimizationCovariance<T> : OptimizationAttribute<NamedMatrix<T>>,
     PolymorphicAttribute<NamedMatrix<T>>(safeTypeOf())
 
+public fun <T> AttributesBuilder<OptimizationProblem<T>>.covariance(covariance: NamedMatrix<T>) {
+    set(OptimizationCovariance(),covariance)
+}
+
 
 public class OptimizationResult<T>() : OptimizationAttribute<Map<Symbol, T>>,
     PolymorphicAttribute<Map<Symbol, T>>(safeTypeOf())
 
+public fun <T> AttributesBuilder<OptimizationProblem<T>>.result(result: Map<Symbol, T>) {
+    set(OptimizationResult(), result)
+}
+
 public val <T> OptimizationProblem<T>.resultOrNull: Map<Symbol, T>? get() = attributes[OptimizationResult()]
 
 public val <T> OptimizationProblem<T>.result: Map<Symbol, T>

kmath-optimization/src/commonMain/kotlin/space/kscience/kmath/optimization/QowOptimizer.kt

@@ -264,10 +264,10 @@ public object QowOptimizer : Optimizer<Double, XYFit> {
             qow = QoWeight(problem, res.freeParameters)
             res = qow.newtonianRun(maxSteps = iterations)
         }
-        val covariance = res.covariance()
+
         return res.problem.withAttributes {
-            set(OptimizationResult(), res.freeParameters)
+            result(res.freeParameters)
-            set(OptimizationCovariance(), covariance)
+            covariance(res.covariance())
         }
     }
 }

kmath-stat/build.gradle.kts

@@ -6,6 +6,7 @@ kscience{
     jvm()
     js()
     native()
+    wasm()
 }
 
 kotlin.sourceSets {

kmath-stat/src/jvmTest/kotlin/space/kscience/kmath/stat/MCScopeTest.kt

@@ -68,7 +68,7 @@ internal class MCScopeTest {
     }
 
 
-    @OptIn(DelicateCoroutinesApi::class)
+    @OptIn(DelicateCoroutinesApi::class, ExperimentalCoroutinesApi::class)
     fun compareResult(test: ATest) {
         val res1 = runBlocking(Dispatchers.Default) { test() }
         val res2 = runBlocking(newSingleThreadContext("test")) { test() }

kmath-symja/src/main/kotlin/space/kscience/kmath/symja/SymjaExpression.kt

@@ -7,8 +7,9 @@ package space.kscience.kmath.symja
 
 import org.matheclipse.core.eval.ExprEvaluator
 import org.matheclipse.core.expression.F
-import space.kscience.kmath.expressions.SpecialDifferentiableExpression
+import space.kscience.attributes.SafeType
 import space.kscience.kmath.expressions.MST
+import space.kscience.kmath.expressions.SpecialDifferentiableExpression
 import space.kscience.kmath.expressions.Symbol
 import space.kscience.kmath.expressions.interpret
 import space.kscience.kmath.operations.NumericAlgebra
@@ -30,6 +31,8 @@ public class SymjaExpression<T : Number, A : NumericAlgebra<T>>(
     public val mst: MST,
     public val evaluator: ExprEvaluator = DEFAULT_EVALUATOR,
 ) : SpecialDifferentiableExpression<T, SymjaExpression<T, A>> {
+    override val type: SafeType<T> get() = algebra.type
+
     override fun invoke(arguments: Map<Symbol, T>): T = mst.interpret(algebra, arguments)
 
     override fun derivativeOrNull(symbols: List<Symbol>): SymjaExpression<T, A> = SymjaExpression(

kmath-tensorflow/src/main/kotlin/space/kscience/kmath/tensorflow/DoubleTensorFlowAlgebra.kt

@@ -10,6 +10,7 @@ import org.tensorflow.Output
 import org.tensorflow.ndarray.NdArray
 import org.tensorflow.op.core.Constant
 import org.tensorflow.types.TFloat64
+import space.kscience.attributes.SafeType
 import space.kscience.kmath.PerformancePitfall
 import space.kscience.kmath.UnstableKMathAPI
 import space.kscience.kmath.expressions.Symbol
@@ -25,9 +26,9 @@ public class DoubleTensorFlowOutput(
     graph: Graph,
     output: Output<TFloat64>,
 ) : TensorFlowOutput<Double, TFloat64>(graph, output) {
+    override val type: SafeType<Double> get() = Float64Field.type
 
     override fun org.tensorflow.Tensor.actualizeTensor(): NdArray<Double> = this as TFloat64
-
 }
 
 internal fun ShapeND.toLongArray(): LongArray = LongArray(size) { get(it).toLong() }

kmath-tensorflow/src/main/kotlin/space/kscience/kmath/tensorflow/IntTensorFlowAlgebra.kt

@@ -10,11 +10,19 @@ import org.tensorflow.Output
 import org.tensorflow.ndarray.NdArray
 import org.tensorflow.types.TInt32
 import org.tensorflow.types.TInt64
+import space.kscience.attributes.SafeType
+import space.kscience.kmath.operations.Int32Ring
+import space.kscience.kmath.operations.Int64Ring
+import space.kscience.kmath.structures.Int32
+import space.kscience.kmath.structures.Int64
 
 public class IntTensorFlowOutput(
     graph: Graph,
     output: Output<TInt32>,
 ) : TensorFlowOutput<Int, TInt32>(graph, output) {
+
+    override val type: SafeType<Int32> get() = Int32Ring.type
+
     override fun org.tensorflow.Tensor.actualizeTensor(): NdArray<Int> = this as TInt32
 }
 
@@ -22,5 +30,7 @@ public class LongTensorFlowOutput(
     graph: Graph,
     output: Output<TInt64>,
 ) : TensorFlowOutput<Long, TInt64>(graph, output) {
+
+    override val type: SafeType<Int64> get() = Int64Ring.type
     override fun org.tensorflow.Tensor.actualizeTensor(): NdArray<Long> = this as TInt64
 }

kmath-tensorflow/src/main/kotlin/space/kscience/kmath/tensorflow/TensorFlowAlgebra.kt

@@ -17,6 +17,7 @@ import org.tensorflow.op.core.*
 import org.tensorflow.types.TInt32
 import org.tensorflow.types.family.TNumber
 import org.tensorflow.types.family.TType
+import space.kscience.attributes.SafeType
 import space.kscience.kmath.PerformancePitfall
 import space.kscience.kmath.UnsafeKMathAPI
 import space.kscience.kmath.UnstableKMathAPI
@@ -39,8 +40,8 @@ public sealed interface TensorFlowTensor<T> : Tensor<T>
 /**
  * Static (eager) in-memory TensorFlow tensor
  */
-@JvmInline
+public class TensorFlowArray<T>(override val type: SafeType<T>, public val tensor: NdArray<T>) : Tensor<T> {
-public value class TensorFlowArray<T>(public val tensor: NdArray<T>) : Tensor<T> {
+
     override val shape: ShapeND get() = ShapeND(tensor.shape().asArray().toIntArray())
 
     @PerformancePitfall
@@ -73,7 +74,7 @@ public abstract class TensorFlowOutput<T, TT : TType>(
 
     internal val actualTensor by lazy {
         Session(graph).use { session ->
-            TensorFlowArray(session.runner().fetch(output).run().first().actualizeTensor())
+            TensorFlowArray(type, session.runner().fetch(output).run().first().actualizeTensor())
         }
     }
 

kmath-tensors/build.gradle.kts

@@ -12,6 +12,7 @@ kscience{
         }
     }
     native()
+    wasm()
 
     dependencies {
         api(projects.kmathCore)

kmath-tensors/src/commonMain/kotlin/space/kscience/kmath/tensors/api/Tensor.kt

@@ -1,10 +0,0 @@
-/*
- * Copyright 2018-2024 KMath contributors.
- * Use of this source code is governed by the Apache 2.0 license that can be found in the license/LICENSE.txt file.
- */
-
-package space.kscience.kmath.tensors.api
-
-import space.kscience.kmath.nd.MutableStructureND
-
-public typealias Tensor<T> = MutableStructureND<T>

kmath-tensors/src/commonMain/kotlin/space/kscience/kmath/tensors/api/TensorAlgebra.kt

@@ -5,11 +5,14 @@
 
 package space.kscience.kmath.tensors.api
 
+import space.kscience.kmath.nd.MutableStructureND
 import space.kscience.kmath.nd.RingOpsND
 import space.kscience.kmath.nd.ShapeND
 import space.kscience.kmath.nd.StructureND
 import space.kscience.kmath.operations.Ring
 
+public typealias Tensor<T> = MutableStructureND<T>
+
 /**
  * Algebra over a ring on [Tensor].
  * For more information: https://proofwiki.org/wiki/Definition:Algebra_over_Ring

kmath-tensors/src/commonMain/kotlin/space/kscience/kmath/tensors/core/LevenbergMarquardtAlgorithm.kt

@@ -6,7 +6,6 @@
 package space.kscience.kmath.tensors.core
 
 import space.kscience.kmath.PerformancePitfall
-import space.kscience.kmath.linear.transposed
 import space.kscience.kmath.nd.*
 import kotlin.math.abs
 import kotlin.math.max
@@ -139,7 +138,7 @@ public fun DoubleTensorAlgebra.levenbergMarquardt(inputData: LMInput): LMResultI
     if (inputData.nargin < 5) {
         weight = fromArray(
             ShapeND(intArrayOf(1, 1)),
-            doubleArrayOf((inputData.realValues.transposed dot inputData.realValues).as1D()[0])
+            doubleArrayOf((inputData.realValues.transposed() dot inputData.realValues).as1D()[0])
         ).as2D()
     }
 
@@ -266,12 +265,12 @@ public fun DoubleTensorAlgebra.levenbergMarquardt(inputData: LMInput): LMResultI
         settings.funcCalls += 1
 
 //        val tmp = deltaY.times(weight)
-        var X2Try = deltaY.as2D().transposed dot deltaY.times(weight) // Chi-squared error criteria
+        var X2Try = deltaY.as2D().transposed() dot deltaY.times(weight) // Chi-squared error criteria
 
         val alpha = 1.0
         if (updateType == 2) { // Quadratic
             // One step of quadratic line update in the h direction for minimum X2
-            val alphaTensor = (jtWdy.transposed dot h) / ((X2Try - x2) / 2.0 + 2 * (jtWdy.transposed dot h))
+            val alphaTensor = (jtWdy.transposed() dot h) / ((X2Try - x2) / 2.0 + 2 * (jtWdy.transposed() dot h))
             h = h dot alphaTensor
             pTry = (p + h).as2D() // update only [idx] elements
             pTry = smallestElementComparison(
@@ -289,7 +288,7 @@ public fun DoubleTensorAlgebra.levenbergMarquardt(inputData: LMInput): LMResultI
             ) // residual error using p_try
             settings.funcCalls += 1
 
-            X2Try = deltaY.as2D().transposed dot deltaY * weight // Chi-squared error criteria
+            X2Try = deltaY.as2D().transposed() dot deltaY * weight // Chi-squared error criteria
         }
 
         val rho = when (updateType) { // Nielsen

kmath-viktor/src/main/kotlin/space/kscience/kmath/viktor/ViktorBuffer.kt

@@ -6,12 +6,16 @@
 package space.kscience.kmath.viktor
 
 import org.jetbrains.bio.viktor.F64FlatArray
+import space.kscience.attributes.SafeType
+import space.kscience.kmath.operations.Float64Field
 import space.kscience.kmath.structures.Buffer
 import space.kscience.kmath.structures.MutableBuffer
 
 @Suppress("NOTHING_TO_INLINE", "OVERRIDE_BY_INLINE")
 @JvmInline
 public value class ViktorBuffer(public val flatArray: F64FlatArray) : MutableBuffer<Double> {
+    override val type: SafeType<Double> get() = Float64Field.type
+
     override val size: Int
         get() = flatArray.length
 

kmath-viktor/src/main/kotlin/space/kscience/kmath/viktor/ViktorStructureND.kt

@@ -6,13 +6,17 @@
 package space.kscience.kmath.viktor
 
 import org.jetbrains.bio.viktor.F64Array
+import space.kscience.attributes.SafeType
 import space.kscience.kmath.PerformancePitfall
 import space.kscience.kmath.nd.ColumnStrides
 import space.kscience.kmath.nd.MutableStructureND
 import space.kscience.kmath.nd.ShapeND
+import space.kscience.kmath.operations.Float64Field
 
 @Suppress("OVERRIDE_BY_INLINE", "NOTHING_TO_INLINE")
 public class ViktorStructureND(public val f64Buffer: F64Array) : MutableStructureND<Double> {
+    override val type: SafeType<Double> get() = Float64Field.type
+
     override val shape: ShapeND get() = ShapeND(f64Buffer.shape)
 
     @OptIn(PerformancePitfall::class)