Optimized performance for Double BufferMatrix product

benchmarks/build.gradle

@@ -1,13 +1,14 @@
 plugins {
     id "java"
-    id "kotlin"
     id "me.champeau.gradle.jmh" version "0.4.7"
+    id 'org.jetbrains.kotlin.jvm'
 }
 
 dependencies {
-    compile project(":kmath-core")
+    api project(":kmath-core")
-    compile project(":kmath-coroutines")
+    api project(":kmath-coroutines")
-    compile project(":kmath-commons")
+    api project(":kmath-commons")
+    implementation "org.jetbrains.kotlin:kotlin-stdlib-jdk8"
     //jmh project(':kmath-core')
 }
 
@@ -16,3 +17,18 @@ jmh{
 }
 
 jmhClasses.dependsOn(compileKotlin)
+repositories {
+    maven { url 'https://dl.bintray.com/kotlin/kotlin-eap' }
+    mavenCentral()
+}
+
+compileKotlin {
+    kotlinOptions {
+        jvmTarget = "1.8"
+    }
+}
+compileTestKotlin {
+    kotlinOptions {
+        jvmTarget = "1.8"
+    }
+}

benchmarks/src/main/kotlin/scientifik/kmath/linear/LinearAlgebraBenchmark.kt

@@ -1,6 +1,5 @@
 package scientifik.kmath.linear
 
-import org.apache.commons.math3.linear.Array2DRowRealMatrix
 import kotlin.random.Random
 import kotlin.system.measureTimeMillis
 
@@ -28,17 +27,13 @@ fun main() {
 
     println("[kmath] Inversion of $n matrices $dim x $dim finished in $inverseTime millis")
 
-    //commons
+    //commons-math
 
     val cmSolver = CMSolver
 
-    val commonsMatrix = Array2DRowRealMatrix(dim, dim)
-    matrix.elements().forEach { (index, value) -> commonsMatrix.setEntry(index[0], index[1], value) }
-
     val commonsTime = measureTimeMillis {
-        val cm = matrix.toCM()
+        val cm = matrix.toCM()             //avoid overhead on conversion
         repeat(n) {
-            //overhead on coversion could be mitigated
             val res = cmSolver.inverse(cm)
         }
     }

benchmarks/src/main/kotlin/scientifik/kmath/linear/MultiplicationBenchmark.kt

@@ -0,0 +1,30 @@
+package scientifik.kmath.linear
+
+import kotlin.random.Random
+import kotlin.system.measureTimeMillis
+
+fun main() {
+    val random = Random(12224)
+    val dim = 1000
+    //creating invertible matrix
+    val matrix1 = Matrix.real(dim, dim) { i, j -> if (i <= j) random.nextDouble() else 0.0 }
+    val matrix2 = Matrix.real(dim, dim) { i, j -> if (i <= j) random.nextDouble() else 0.0 }
+
+//    //warmup
+//    matrix1 dot matrix2
+
+    val cmMatrix1 = matrix1.toCM()
+    val cmMatrix2 = matrix2.toCM()
+
+    val cmTime = measureTimeMillis {
+        cmMatrix1 dot cmMatrix2
+    }
+
+    println("CM implementation time: $cmTime")
+
+    val genericTime = measureTimeMillis {
+        val res = matrix1 dot matrix2
+    }
+
+    println("Generic implementation time: $genericTime")
+}

build.gradle.kts

@@ -28,7 +28,7 @@ allprojects {
     apply(plugin = "com.jfrog.artifactory")
 
     group = "scientifik"
-    version = "0.0.3-dev-2"
+    version = "0.0.3-dev-3"
 
     repositories {
         maven("https://dl.bintray.com/kotlin/kotlin-eap")

kmath-commons/build.gradle.kts

@@ -2,14 +2,11 @@ plugins {
     kotlin("jvm")
 }
 
+description = "Commons math binding for kmath"
+
 dependencies {
     api(project(":kmath-core"))
     api("org.apache.commons:commons-math3:3.6.1")
     testImplementation("org.jetbrains.kotlin:kotlin-test")
     testImplementation("org.jetbrains.kotlin:kotlin-test-junit")
 }
-
-//dependencies {
-////    compile(project(":kmath-core"))
-////    //compile project(":kmath-coroutines")
-////}

kmath-commons/src/main/kotlin/scientifik/kmath/linear/CMMatrix.kt

@@ -55,7 +55,7 @@ object CMMatrixContext : MatrixContext<Double, RealField> {
         return CMVector(ArrayRealVector(array))
     }
 
-    override fun Matrix<Double>.dot(other: Matrix<Double>): Matrix<Double> = this.toCM().dot(other.toCM())
+    override fun Matrix<Double>.dot(other: Matrix<Double>): Matrix<Double> = CMMatrix(this.toCM().origin.multiply(other.toCM().origin))
 }
 
 operator fun CMMatrix.plus(other: CMMatrix): CMMatrix = CMMatrix(this.origin.add(other.origin))

kmath-core/build.gradle

@@ -3,14 +3,11 @@ plugins {
 }
 
 kotlin {
-    targets {
+    jvm {
-        fromPreset(presets.jvm, 'jvm')
+        compilations["main"].kotlinOptions.jvmTarget = "1.8"
-        fromPreset(presets.js, 'js')
-        // For ARM, preset should be changed to presets.iosArm32 or presets.iosArm64
-        // For Linux, preset should be changed to e.g. presets.linuxX64
-        // For MacOS, preset should be changed to e.g. presets.macosX64
-        //fromPreset(presets.mingwX64, 'mingw')
     }
+    js()
+
     sourceSets {
         commonMain {
             dependencies {

kmath-core/src/commonMain/kotlin/scientifik/kmath/linear/BufferMatrix.kt

@@ -0,0 +1,97 @@
+package scientifik.kmath.linear
+
+import scientifik.kmath.operations.Ring
+import scientifik.kmath.structures.*
+
+/**
+ * Basic implementation of Matrix space based on [NDStructure]
+ */
+class BufferMatrixContext<T : Any, R : Ring<T>>(
+    override val elementContext: R,
+    private val bufferFactory: BufferFactory<T>
+) : MatrixContext<T, R> {
+
+    override fun produce(rows: Int, columns: Int, initializer: (i: Int, j: Int) -> T): BufferMatrix<T> {
+        val buffer = bufferFactory(rows * columns) { offset -> initializer(offset / columns, offset % columns) }
+        return BufferMatrix(rows, columns, buffer)
+    }
+
+    override fun point(size: Int, initializer: (Int) -> T): Point<T> = bufferFactory(size, initializer)
+}
+
+class BufferMatrix<T : Any>(
+    override val rowNum: Int,
+    override val colNum: Int,
+    val buffer: Buffer<out T>,
+    override val features: Set<MatrixFeature> = emptySet()
+) : Matrix<T> {
+
+    init {
+        if (buffer.size != rowNum * colNum) {
+            error("Dimension mismatch for matrix structure")
+        }
+    }
+
+
+    override val shape: IntArray get() = intArrayOf(rowNum, colNum)
+
+    override fun get(index: IntArray): T = get(index[0], index[1])
+
+    override fun get(i: Int, j: Int): T = buffer[i * colNum + j]
+
+    override fun elements(): Sequence<Pair<IntArray, T>> = sequence {
+        for (i in 0 until rowNum) {
+            for (j in 0 until colNum) {
+                yield(intArrayOf(i, j) to get(i, j))
+            }
+        }
+    }
+
+    override fun equals(other: Any?): Boolean {
+        if (this === other) return true
+        return when (other) {
+            is NDStructure<*> -> return NDStructure.equals(this, other)
+            else -> false
+        }
+    }
+
+    override fun hashCode(): Int {
+        var result = buffer.hashCode()
+        result = 31 * result + features.hashCode()
+        return result
+    }
+
+    override fun toString(): String {
+        return if (rowNum <= 5 && colNum <= 5) {
+            "Matrix(rowsNum = $rowNum, colNum = $colNum, features=$features)\n" +
+                    rows.asSequence().joinToString(prefix = "(", postfix = ")", separator = "\n ") {
+                        it.asSequence().joinToString(separator = "\t") { it.toString() }
+                    }
+        } else {
+            "Matrix(rowsNum = $rowNum, colNum = $colNum, features=$features)"
+        }
+    }
+}
+
+/**
+ * Optimized dot product for real matrices
+ */
+infix fun BufferMatrix<Double>.dot(other: BufferMatrix<Double>): BufferMatrix<Double> {
+    if (this.colNum != other.rowNum) error("Matrix dot operation dimension mismatch: ($rowNum, $colNum) x (${other.rowNum}, ${other.colNum})")
+
+    val array = DoubleArray(this.rowNum * other.colNum)
+
+    val a = this.buffer.array
+    val b = other.buffer.array
+
+    for (i in (0 until rowNum)) {
+        for (j in (0 until other.colNum)) {
+            for (k in (0 until colNum)) {
+                array[i * other.colNum + j] += a[i * colNum + k] * b[k * other.colNum + j]
+            }
+        }
+    }
+
+    val buffer = DoubleBuffer(array)
+    return BufferMatrix(rowNum, other.colNum, buffer)
+}

kmath-core/src/commonMain/kotlin/scientifik/kmath/linear/Matrix.kt

@@ -75,13 +75,13 @@ interface MatrixContext<T : Any, R : Ring<T>> {
         /**
          * Non-boxing double matrix
          */
-        val real: MatrixContext<Double, RealField> = StructureMatrixContext(RealField, DoubleBufferFactory)
+        val real = BufferMatrixContext(RealField, DoubleBufferFactory)
 
         /**
          * A structured matrix with custom buffer
          */
         fun <T : Any, R : Ring<T>> buffered(ring: R, bufferFactory: BufferFactory<T> = ::boxing): MatrixContext<T, R> =
-            StructureMatrixContext(ring, bufferFactory)
+            BufferMatrixContext(ring, bufferFactory)
 
         /**
          * Automatic buffered matrix, unboxed if it is possible
@@ -152,11 +152,10 @@ interface Matrix<T : Any> : NDStructure<T> {
          * Build a square matrix from given elements.
          */
         fun <T : Any> build(vararg elements: T): Matrix<T> {
-            val buffer = elements.asBuffer()
             val size: Int = sqrt(elements.size.toDouble()).toInt()
             if (size * size != elements.size) error("The number of elements ${elements.size} is not a full square")
-            val structure = Mutable2DStructure(size, size, buffer)
+            val buffer = elements.asBuffer()
-            return StructureMatrix(structure)
+            return BufferMatrix(size, size, buffer)
         }
     }
 }

kmath-core/src/commonMain/kotlin/scientifik/kmath/linear/StructureMatrix.kt

@@ -1,74 +0,0 @@
-package scientifik.kmath.linear
-
-import scientifik.kmath.operations.Ring
-import scientifik.kmath.structures.*
-
-/**
- * Basic implementation of Matrix space based on [NDStructure]
- */
-class StructureMatrixContext<T : Any, R : Ring<T>>(
-    override val elementContext: R,
-    private val bufferFactory: BufferFactory<T>
-) : MatrixContext<T, R> {
-
-    override fun produce(rows: Int, columns: Int, initializer: (i: Int, j: Int) -> T): Matrix<T> {
-        val structure =
-            ndStructure(intArrayOf(rows, columns), bufferFactory) { index -> initializer(index[0], index[1]) }
-        return StructureMatrix(structure)
-    }
-
-    override fun point(size: Int, initializer: (Int) -> T): Point<T> = bufferFactory(size, initializer)
-}
-
-class StructureMatrix<T : Any>(
-    val structure: NDStructure<out T>,
-    override val features: Set<MatrixFeature> = emptySet()
-) : Matrix<T> {
-
-    init {
-        if (structure.shape.size != 2) {
-            error("Dimension mismatch for matrix structure")
-        }
-    }
-
-    override val rowNum: Int
-        get() = structure.shape[0]
-    override val colNum: Int
-        get() = structure.shape[1]
-
-
-    override val shape: IntArray get() = structure.shape
-
-    override fun get(index: IntArray): T = structure[index]
-
-    override fun get(i: Int, j: Int): T = structure[i, j]
-
-    override fun elements(): Sequence<Pair<IntArray, T>> = structure.elements()
-
-    override fun equals(other: Any?): Boolean {
-        if (this === other) return true
-        return when (other) {
-            is NDStructure<*> -> return NDStructure.equals(this, other)
-            else -> false
-        }
-    }
-
-    override fun hashCode(): Int {
-        var result = structure.hashCode()
-        result = 31 * result + features.hashCode()
-        return result
-    }
-
-    override fun toString(): String {
-        return if (rowNum <= 5 && colNum <= 5) {
-            "Matrix(rowsNum = $rowNum, colNum = $colNum, features=$features)\n" +
-                    rows.asSequence().joinToString(prefix = "(", postfix = ")", separator = "\n ") {
-                        it.asSequence().joinToString(separator = "\t") { it.toString() }
-                    }
-        } else {
-            "Matrix(rowsNum = $rowNum, colNum = $colNum, features=$features)"
-        }
-    }
-
-
-}

kmath-core/src/commonMain/kotlin/scientifik/kmath/structures/Buffers.kt

@@ -97,7 +97,7 @@ interface MutableBuffer<T> : Buffer<T> {
 }
 
 
-inline class ListBuffer<T>(private val list: List<T>) : Buffer<T> {
+inline class ListBuffer<T>(val list: List<T>) : Buffer<T> {
 
     override val size: Int
         get() = list.size
@@ -109,7 +109,7 @@ inline class ListBuffer<T>(private val list: List<T>) : Buffer<T> {
 
 fun <T> List<T>.asBuffer() = ListBuffer(this)
 
-inline class MutableListBuffer<T>(private val list: MutableList<T>) : MutableBuffer<T> {
+inline class MutableListBuffer<T>(val list: MutableList<T>) : MutableBuffer<T> {
 
     override val size: Int
         get() = list.size
@@ -142,7 +142,7 @@ class ArrayBuffer<T>(private val array: Array<T>) : MutableBuffer<T> {
 
 fun <T> Array<T>.asBuffer() = ArrayBuffer(this)
 
-inline class DoubleBuffer(private val array: DoubleArray) : MutableBuffer<Double> {
+inline class DoubleBuffer(val array: DoubleArray) : MutableBuffer<Double> {
     override val size: Int get() = array.size
 
     override fun get(index: Int): Double = array[index]
@@ -157,9 +157,19 @@ inline class DoubleBuffer(private val array: DoubleArray) : MutableBuffer<Double
 
 }
 
+/**
+ * Transform buffer of doubles into array for high performance operations
+ */
+val Buffer<out Double>.array: DoubleArray
+    get() = if (this is DoubleBuffer) {
+        array
+    } else {
+        DoubleArray(size) { get(it) }
+    }
+
 fun DoubleArray.asBuffer() = DoubleBuffer(this)
 
-inline class ShortBuffer(private val array: ShortArray) : MutableBuffer<Short> {
+inline class ShortBuffer(val array: ShortArray) : MutableBuffer<Short> {
     override val size: Int get() = array.size
 
     override fun get(index: Int): Short = array[index]
@@ -176,7 +186,7 @@ inline class ShortBuffer(private val array: ShortArray) : MutableBuffer<Short> {
 
 fun ShortArray.asBuffer() = ShortBuffer(this)
 
-inline class IntBuffer(private val array: IntArray) : MutableBuffer<Int> {
+inline class IntBuffer(val array: IntArray) : MutableBuffer<Int> {
     override val size: Int get() = array.size
 
     override fun get(index: Int): Int = array[index]
@@ -193,7 +203,7 @@ inline class IntBuffer(private val array: IntArray) : MutableBuffer<Int> {
 
 fun IntArray.asBuffer() = IntBuffer(this)
 
-inline class LongBuffer(private val array: LongArray) : MutableBuffer<Long> {
+inline class LongBuffer(val array: LongArray) : MutableBuffer<Long> {
     override val size: Int get() = array.size
 
     override fun get(index: Int): Long = array[index]
@@ -210,7 +220,7 @@ inline class LongBuffer(private val array: LongArray) : MutableBuffer<Long> {
 
 fun LongArray.asBuffer() = LongBuffer(this)
 
-inline class ReadOnlyBuffer<T>(private val buffer: MutableBuffer<T>) : Buffer<T> {
+inline class ReadOnlyBuffer<T>(val buffer: MutableBuffer<T>) : Buffer<T> {
     override val size: Int get() = buffer.size
 
     override fun get(index: Int): T = buffer.get(index)

kmath-core/src/commonTest/kotlin/scientifik/kmath/linear/MatrixTest.kt

@@ -24,7 +24,7 @@ class MatrixTest {
     fun testTranspose() {
         val matrix = MatrixContext.real.one(3, 3)
         val transposed = matrix.transpose()
-        assertEquals((matrix as StructureMatrix).structure, (transposed as StructureMatrix).structure)
+        assertEquals((matrix as BufferMatrix).buffer, (transposed as BufferMatrix).buffer)
         assertEquals(matrix, transposed)
     }
 

kmath-koma/build.gradle.kts

@@ -0,0 +1,53 @@
+plugins {
+    id("kotlin-multiplatform")
+}
+
+repositories {
+    maven("http://dl.bintray.com/kyonifer/maven")
+}
+
+kotlin {
+    jvm {
+        compilations["main"].kotlinOptions.jvmTarget = "1.8"
+    }
+    js()
+
+    sourceSets {
+
+        val commonMain by getting {
+            dependencies {
+                api(project(":kmath-core"))
+                implementation("com.kyonifer:koma-core-api-common:0.12")
+                implementation("org.jetbrains.kotlin:kotlin-stdlib-common")
+            }
+        }
+        val commonTest by getting {
+            dependencies {
+                implementation(kotlin("test-common"))
+                implementation(kotlin("test-annotations-common"))
+            }
+        }
+        val jvmMain by getting {
+            dependencies {
+                implementation(kotlin("stdlib-jdk8"))
+            }
+        }
+        val jvmTest by getting {
+            dependencies {
+                implementation(kotlin("test"))
+                implementation(kotlin("test-junit"))
+                implementation("com.kyonifer:koma-core-ejml:0.12")
+            }
+        }
+        val jsMain by getting {
+            dependencies {
+                implementation(kotlin("stdlib-js"))
+            }
+        }
+        val jsTest by getting {
+            dependencies {
+                implementation(kotlin("test-js"))
+            }
+        }
+    }
+}

kmath-koma/src/commonMain/kotlin/scientifik.kmath.linear/KomaMatrix.kt

@@ -0,0 +1,27 @@
+package scientifik.kmath.linear
+
+import scientifik.kmath.operations.Ring
+
+class KomaMatrixContext<T: Any, R: Ring<T>> : MatrixContext<T,R> {
+    override val elementContext: R
+        get() = TODO("not implemented") //To change initializer of created properties use File | Settings | File Templates.
+
+    override fun produce(rows: Int, columns: Int, initializer: (i: Int, j: Int) -> T): Matrix<T> {
+        TODO("not implemented") //To change body of created functions use File | Settings | File Templates.
+    }
+
+    override fun point(size: Int, initializer: (Int) -> T): Point<T> {
+        TODO("not implemented") //To change body of created functions use File | Settings | File Templates.
+    }
+
+}
+
+inline class KomaMatrix<T : Any>(val matrix: koma.matrix.Matrix<T>) : Matrix<T> {
+    override val rowNum: Int get() = matrix.numRows()
+    override val colNum: Int get() = matrix.numCols()
+    override val features: Set<MatrixFeature> get() = emptySet()
+
+    @Suppress("OVERRIDE_BY_INLINE")
+    override inline fun get(i: Int, j: Int): T = matrix.getGeneric(i, j)
+
+}

settings.gradle.kts

@@ -2,6 +2,8 @@ pluginManagement {
     repositories {
         mavenCentral()
         maven("https://plugins.gradle.org/m2/")
+        maven { setUrl("https://dl.bintray.com/kotlin/kotlin-eap") }
+        maven { setUrl("https://plugins.gradle.org/m2/") }
     }
 }
 
@@ -13,5 +15,6 @@ include(
     ":kmath-io",
     ":kmath-coroutines",
     ":kmath-commons",
+    ":kmath-koma",
     ":benchmarks"
 )