Multik integration finished (for now)

kmath-core/src/commonMain/kotlin/space/kscience/kmath/operations/DoubleBufferOps.kt

@@ -17,12 +17,15 @@ import kotlin.math.*
 /**
  * [ExtendedFieldOps] over [DoubleBuffer].
  */
-public abstract class DoubleBufferOps :
+public abstract class DoubleBufferOps : BufferAlgebra<Double, DoubleField>, ExtendedFieldOps<Buffer<Double>>,
-    BufferAlgebra<Double, DoubleField>, ExtendedFieldOps<Buffer<Double>>, Norm<Buffer<Double>, Double> {
+    Norm<Buffer<Double>, Double> {
 
     override val elementAlgebra: DoubleField get() = DoubleField
     override val bufferFactory: BufferFactory<Double> get() = ::DoubleBuffer
 
+    override fun Buffer<Double>.map(block: DoubleField.(Double) -> Double): DoubleBuffer =
+        mapInline { DoubleField.block(it) }
+
     @UnstableKMathAPI
     override fun unaryOperationFunction(operation: String): (arg: Buffer<Double>) -> Buffer<Double> =
         super<ExtendedFieldOps>.unaryOperationFunction(operation)
@@ -87,8 +90,7 @@ public abstract class DoubleBufferOps :
             val aArray = left.array
             val bArray = right.array
             DoubleBuffer(DoubleArray(left.size) { aArray[it] * bArray[it] })
-        } else
+        } else DoubleBuffer(DoubleArray(left.size) { left[it] * right[it] })
-            DoubleBuffer(DoubleArray(left.size) { left[it] * right[it] })
     }
 
     override fun divide(left: Buffer<Double>, right: Buffer<Double>): DoubleBuffer {

kmath-multik/src/main/kotlin/space/kscience/kmath/multik/MultikDoubleAlgebra.kt

@@ -1,134 +1,47 @@
 package space.kscience.kmath.multik
 
-import org.jetbrains.kotlinx.multik.ndarray.data.DN
 import org.jetbrains.kotlinx.multik.ndarray.data.DataType
 import space.kscience.kmath.nd.StructureND
 import space.kscience.kmath.operations.DoubleField
-import space.kscience.kmath.tensors.api.AnalyticTensorAlgebra
+import space.kscience.kmath.operations.ExponentialOperations
-import space.kscience.kmath.tensors.api.LinearOpsTensorAlgebra
+import space.kscience.kmath.operations.TrigonometricOperations
-import space.kscience.kmath.tensors.api.Tensor
 
 public object MultikDoubleAlgebra : MultikDivisionTensorAlgebra<Double, DoubleField>(),
-    AnalyticTensorAlgebra<Double, DoubleField>, LinearOpsTensorAlgebra<Double, DoubleField> {
+    TrigonometricOperations<StructureND<Double>>, ExponentialOperations<StructureND<Double>> {
     override val elementAlgebra: DoubleField get() = DoubleField
     override val type: DataType get() = DataType.DoubleDataType
 
-    override fun StructureND<Double>.mean(): Double = multikStat.mean(asMultik().array)
+    override fun sin(arg: StructureND<Double>): MultikTensor<Double> = multikMath.mathEx.sin(arg.asMultik().array).wrap()
 
-    override fun StructureND<Double>.mean(dim: Int, keepDim: Boolean): Tensor<Double> =
+    override fun cos(arg: StructureND<Double>): MultikTensor<Double> = multikMath.mathEx.cos(arg.asMultik().array).wrap()
-        multikStat.mean<Double,DN, DN>(asMultik().array, dim).wrap()
 
-    override fun StructureND<Double>.std(): Double {
+    override fun tan(arg: StructureND<Double>): MultikTensor<Double> = sin(arg) / cos(arg)
-        TODO("Not yet implemented")
+
+    override fun asin(arg: StructureND<Double>): MultikTensor<Double>  = arg.map { asin(it) }
+
+    override fun acos(arg: StructureND<Double>): MultikTensor<Double> = arg.map { acos(it) }
+
+    override fun atan(arg: StructureND<Double>): MultikTensor<Double> = arg.map { atan(it) }
+
+    override fun exp(arg: StructureND<Double>): MultikTensor<Double> = multikMath.mathEx.exp(arg.asMultik().array).wrap()
+
+    override fun ln(arg: StructureND<Double>): MultikTensor<Double> = multikMath.mathEx.log(arg.asMultik().array).wrap()
+
+    override fun sinh(arg: StructureND<Double>): MultikTensor<Double> = (exp(arg) - exp(-arg)) / 2.0
+
+    override fun cosh(arg: StructureND<Double>): MultikTensor<Double> = (exp(arg) + exp(-arg)) / 2.0
+
+    override fun tanh(arg: StructureND<Double>): MultikTensor<Double> {
+        val expPlus = exp(arg)
+        val expMinus = exp(-arg)
+        return (expPlus - expMinus) / (expPlus + expMinus)
     }
 
-    override fun StructureND<Double>.std(dim: Int, keepDim: Boolean): Tensor<Double> {
+    override fun asinh(arg: StructureND<Double>): MultikTensor<Double> = arg.map { asinh(it) }
-        TODO("Not yet implemented")
-    }
 
-    override fun StructureND<Double>.variance(): Double {
+    override fun acosh(arg: StructureND<Double>): MultikTensor<Double> = arg.map { acosh(it) }
-        TODO("Not yet implemented")
-    }
 
-    override fun StructureND<Double>.variance(dim: Int, keepDim: Boolean): Tensor<Double> {
+    override fun atanh(arg: StructureND<Double>): MultikTensor<Double> = arg.map { atanh(it) }
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.exp(): Tensor<Double> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.ln(): Tensor<Double> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.sqrt(): Tensor<Double> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.cos(): Tensor<Double> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.acos(): Tensor<Double> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.cosh(): Tensor<Double> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.acosh(): Tensor<Double> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.sin(): Tensor<Double> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.asin(): Tensor<Double> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.sinh(): Tensor<Double> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.asinh(): Tensor<Double> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.tan(): Tensor<Double> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.atan(): Tensor<Double> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.tanh(): Tensor<Double> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.atanh(): Tensor<Double> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.ceil(): Tensor<Double> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.floor(): Tensor<Double> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.det(): Tensor<Double> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.inv(): Tensor<Double> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.cholesky(): Tensor<Double> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.qr(): Pair<Tensor<Double>, Tensor<Double>> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.lu(): Triple<Tensor<Double>, Tensor<Double>, Tensor<Double>> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.svd(): Triple<Tensor<Double>, Tensor<Double>, Tensor<Double>> {
-        TODO("Not yet implemented")
-    }
-
-    override fun StructureND<Double>.symEig(): Pair<Tensor<Double>, Tensor<Double>> {
-        TODO("Not yet implemented")
-    }
 }
 
 public val Double.Companion.multikAlgebra: MultikTensorAlgebra<Double, DoubleField> get() = MultikDoubleAlgebra

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

@@ -75,9 +75,9 @@ public open class DoubleTensorAlgebra :
     override fun zip(
         left: StructureND<Double>,
         right: StructureND<Double>,
-        transform: DoubleField.(Double, Double) -> Double
+        transform: DoubleField.(Double, Double) -> Double,
     ): DoubleTensor {
-        require(left.shape.contentEquals(right.shape)){
+        require(left.shape.contentEquals(right.shape)) {
             "The shapes in zip are not equal: left - ${left.shape}, right - ${right.shape}"
         }
         val leftTensor = left.tensor
@@ -422,14 +422,11 @@ public open class DoubleTensorAlgebra :
 
         for ((res, ab) in resTensor.matrixSequence().zip(newThis.matrixSequence().zip(newOther.matrixSequence()))) {
             val (a, b) = ab
-            dotHelper(a.as2D(), b.as2D(), res.as2D(), l, m1, n)
+            dotTo(a.as2D(), b.as2D(), res.as2D(), l, m1, n)
         }
 
         if (penultimateDim) {
-            return resTensor.view(
+            return resTensor.view(resTensor.shape.dropLast(2).toIntArray() + intArrayOf(resTensor.shape.last()))
-                resTensor.shape.dropLast(2).toIntArray() +
-                        intArrayOf(resTensor.shape.last())
-            )
         }
         if (lastDim) {
             return resTensor.view(resTensor.shape.dropLast(1).toIntArray())
@@ -441,7 +438,7 @@ public open class DoubleTensorAlgebra :
         diagonalEntries: Tensor<Double>,
         offset: Int,
         dim1: Int,
-        dim2: Int
+        dim2: Int,
     ): DoubleTensor {
         val n = diagonalEntries.shape.size
         val d1 = minusIndexFrom(n + 1, dim1)
@@ -577,13 +574,13 @@ public open class DoubleTensorAlgebra :
      */
     public fun Tensor<Double>.rowsByIndices(indices: IntArray): DoubleTensor = stack(indices.map { this[it] })
 
-    internal inline fun StructureND<Double>.fold(foldFunction: (DoubleArray) -> Double): Double =
+    private inline fun StructureND<Double>.fold(foldFunction: (DoubleArray) -> Double): Double =
         foldFunction(tensor.copyArray())
 
-    internal inline fun <reified R: Any> StructureND<Double>.foldDim(
+    private inline fun <reified R : Any> StructureND<Double>.foldDim(
-        foldFunction: (DoubleArray) -> R,
         dim: Int,
         keepDim: Boolean,
+        foldFunction: (DoubleArray) -> R,
     ): BufferedTensor<R> {
         check(dim < dimension) { "Dimension $dim out of range $dimension" }
         val resShape = if (keepDim) {
@@ -592,7 +589,7 @@ public open class DoubleTensorAlgebra :
             shape.take(dim).toIntArray() + shape.takeLast(dimension - dim - 1).toIntArray()
         }
         val resNumElements = resShape.reduce(Int::times)
-        val init = foldFunction(DoubleArray(1){0.0})
+        val init = foldFunction(DoubleArray(1) { 0.0 })
         val resTensor = BufferedTensor(resShape,
             MutableBuffer.auto(resNumElements) { init }, 0)
         for (index in resTensor.indices) {
@@ -608,66 +605,59 @@ public open class DoubleTensorAlgebra :
     override fun StructureND<Double>.sum(): Double = tensor.fold { it.sum() }
 
     override fun StructureND<Double>.sum(dim: Int, keepDim: Boolean): DoubleTensor =
-        foldDim({ x -> x.sum() }, dim, keepDim).toDoubleTensor()
+        foldDim(dim, keepDim) { x -> x.sum() }.toDoubleTensor()
 
     override fun StructureND<Double>.min(): Double = this.fold { it.minOrNull()!! }
 
     override fun StructureND<Double>.min(dim: Int, keepDim: Boolean): DoubleTensor =
-        foldDim({ x -> x.minOrNull()!! }, dim, keepDim).toDoubleTensor()
+        foldDim(dim, keepDim) { x -> x.minOrNull()!! }.toDoubleTensor()
 
     override fun StructureND<Double>.max(): Double = this.fold { it.maxOrNull()!! }
 
     override fun StructureND<Double>.max(dim: Int, keepDim: Boolean): DoubleTensor =
-        foldDim({ x -> x.maxOrNull()!! }, dim, keepDim).toDoubleTensor()
+        foldDim(dim, keepDim) { x -> x.maxOrNull()!! }.toDoubleTensor()
 
 
     override fun StructureND<Double>.argMax(dim: Int, keepDim: Boolean): IntTensor =
-        foldDim({ x ->
+        foldDim(dim, keepDim) { x ->
             x.withIndex().maxByOrNull { it.value }?.index!!
-        }, dim, keepDim).toIntTensor()
+        }.toIntTensor()
 
 
     override fun StructureND<Double>.mean(): Double = this.fold { it.sum() / tensor.numElements }
 
-    override fun StructureND<Double>.mean(dim: Int, keepDim: Boolean): DoubleTensor =
+    override fun StructureND<Double>.mean(dim: Int, keepDim: Boolean): DoubleTensor = foldDim(dim, keepDim) { arr ->
-        foldDim(
+        check(dim < dimension) { "Dimension $dim out of range $dimension" }
-            { arr ->
+        arr.sum() / shape[dim]
-                check(dim < dimension) { "Dimension $dim out of range $dimension" }
+    }.toDoubleTensor()
-                arr.sum() / shape[dim]
-            },
-            dim,
-            keepDim
-        ).toDoubleTensor()
 
-    override fun StructureND<Double>.std(): Double = this.fold { arr ->
+    override fun StructureND<Double>.std(): Double = fold { arr ->
         val mean = arr.sum() / tensor.numElements
         sqrt(arr.sumOf { (it - mean) * (it - mean) } / (tensor.numElements - 1))
     }
 
     override fun StructureND<Double>.std(dim: Int, keepDim: Boolean): DoubleTensor = foldDim(
-        { arr ->
-            check(dim < dimension) { "Dimension $dim out of range $dimension" }
-            val mean = arr.sum() / shape[dim]
-            sqrt(arr.sumOf { (it - mean) * (it - mean) } / (shape[dim] - 1))
-        },
         dim,
         keepDim
-    ).toDoubleTensor()
+    ) { arr ->
+        check(dim < dimension) { "Dimension $dim out of range $dimension" }
+        val mean = arr.sum() / shape[dim]
+        sqrt(arr.sumOf { (it - mean) * (it - mean) } / (shape[dim] - 1))
+    }.toDoubleTensor()
 
-    override fun StructureND<Double>.variance(): Double = this.fold { arr ->
+    override fun StructureND<Double>.variance(): Double = fold { arr ->
         val mean = arr.sum() / tensor.numElements
         arr.sumOf { (it - mean) * (it - mean) } / (tensor.numElements - 1)
     }
 
     override fun StructureND<Double>.variance(dim: Int, keepDim: Boolean): DoubleTensor = foldDim(
-        { arr ->
-            check(dim < dimension) { "Dimension $dim out of range $dimension" }
-            val mean = arr.sum() / shape[dim]
-            arr.sumOf { (it - mean) * (it - mean) } / (shape[dim] - 1)
-        },
         dim,
         keepDim
-    ).toDoubleTensor()
+    ) { arr ->
+        check(dim < dimension) { "Dimension $dim out of range $dimension" }
+        val mean = arr.sum() / shape[dim]
+        arr.sumOf { (it - mean) * (it - mean) } / (shape[dim] - 1)
+    }.toDoubleTensor()
 
     private fun cov(x: DoubleTensor, y: DoubleTensor): Double {
         val n = x.shape[0]
@@ -699,19 +689,14 @@ public open class DoubleTensorAlgebra :
         return resTensor
     }
 
-    @OptIn(PerformancePitfall::class)
     override fun StructureND<Double>.exp(): DoubleTensor = tensor.map { exp(it) }
 
-    @OptIn(PerformancePitfall::class)
     override fun StructureND<Double>.ln(): DoubleTensor = tensor.map { ln(it) }
 
-    @OptIn(PerformancePitfall::class)
     override fun StructureND<Double>.sqrt(): DoubleTensor = tensor.map { sqrt(it) }
 
-    @OptIn(PerformancePitfall::class)
     override fun StructureND<Double>.cos(): DoubleTensor = tensor.map { cos(it) }
 
-    @OptIn(PerformancePitfall::class)
     override fun StructureND<Double>.acos(): DoubleTensor = tensor.map { acos(it) }
 
     override fun StructureND<Double>.cosh(): DoubleTensor = tensor.map { cosh(it) }

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

@@ -53,7 +53,7 @@ internal val <T> BufferedTensor<T>.matrices: VirtualBuffer<BufferedTensor<T>>
 
 internal fun <T> BufferedTensor<T>.matrixSequence(): Sequence<BufferedTensor<T>> = matrices.asSequence()
 
-internal fun dotHelper(
+internal fun dotTo(
     a: MutableStructure2D<Double>,
     b: MutableStructure2D<Double>,
     res: MutableStructure2D<Double>,