Linear interpolation

README.md

@@ -1,9 +1,12 @@
+[![JetBrains Research](https://jb.gg/badges/research.svg)](https://confluence.jetbrains.com/display/ALL/JetBrains+on+GitHub)
+[![DOI](https://zenodo.org/badge/129486382.svg)](https://zenodo.org/badge/latestdoi/129486382)
+
+![Gradle build](https://github.com/mipt-npm/kmath/workflows/Gradle%20build/badge.svg)
+
 Bintray:        [ ![Download](https://api.bintray.com/packages/mipt-npm/scientifik/kmath-core/images/download.svg) ](https://bintray.com/mipt-npm/scientifik/kmath-core/_latestVersion)
 
 Bintray-dev:    [ ![Download](https://api.bintray.com/packages/mipt-npm/dev/kmath-core/images/download.svg) ](https://bintray.com/mipt-npm/scientifik/kmath-core/_latestVersion)
 
-[![DOI](https://zenodo.org/badge/129486382.svg)](https://zenodo.org/badge/latestdoi/129486382)
-
 # KMath
 Could be pronounced as `key-math`.
 The Kotlin MATHematics library is intended as a Kotlin-based analog to Python's `numpy` library. In contrast to `numpy` and `scipy` it is modular and has a lightweight core.

build.gradle.kts

@@ -1,5 +1,5 @@
 plugins {
-    id("scientifik.publish") version "0.2.6" apply false
+    id("scientifik.publish") version "0.3.1" apply false
 }
 
 val kmathVersion by extra("0.1.4-dev-1")

kmath-core/src/commonMain/kotlin/scientifik/kmath/operations/Algebra.kt

@@ -6,14 +6,14 @@ annotation class KMathContext
 /**
  * Marker interface for any algebra
  */
-interface Algebra
+interface Algebra<T>
 
-inline operator fun <T : Algebra, R> T.invoke(block: T.() -> R): R = run(block)
+inline operator fun <T : Algebra<*>, R> T.invoke(block: T.() -> R): R = run(block)
 
 /**
  * Space-like operations without neutral element
  */
-interface SpaceOperations<T> : Algebra {
+interface SpaceOperations<T> : Algebra<T> {
     /**
      * Addition operation for two context elements
      */

kmath-core/src/commonMain/kotlin/scientifik/kmath/operations/AlgebraExtensions.kt

@@ -1,4 +1,13 @@
 package scientifik.kmath.operations
 
-fun <T> Space<T>.sum(data : Iterable<T>): T = data.fold(zero) { left, right -> add(left,right) }
-fun <T> Space<T>.sum(data : Sequence<T>): T = data.fold(zero) { left, right -> add(left, right) }
+fun <T> Space<T>.sum(data: Iterable<T>): T = data.fold(zero) { left, right -> add(left, right) }
+fun <T> Space<T>.sum(data: Sequence<T>): T = data.fold(zero) { left, right -> add(left, right) }
+
+//TODO optimized power operation
+fun <T> RingOperations<T>.power(arg: T, power: Int): T {
+    var res = arg
+    repeat(power - 1) {
+        res *= arg
+    }
+    return res
+}

kmath-core/src/commonMain/kotlin/scientifik/kmath/operations/OptionalOperations.kt

@@ -29,7 +29,7 @@ fun <T : MathElement<out TrigonometricOperations<T>>> ctg(arg: T): T = arg.conte
 /**
  * A context extension to include power operations like square roots, etc
  */
-interface PowerOperations<T> {
+interface PowerOperations<T> : Algebra<T> {
     fun power(arg: T, pow: Number): T
     fun sqrt(arg: T) = power(arg, 0.5)
 
@@ -42,7 +42,7 @@ fun <T : MathElement<out PowerOperations<T>>> sqr(arg: T): T = arg pow 2.0
 
 /* Exponential */
 
-interface ExponentialOperations<T> {
+interface ExponentialOperations<T>: Algebra<T> {
     fun exp(arg: T): T
     fun ln(arg: T): T
 }

kmath-functions/build.gradle.kts

@@ -1,23 +1,11 @@
 plugins {
     id("scientifik.mpp")
-    //id("scientifik.atomic")
 }
 
 kotlin.sourceSets {
     commonMain {
         dependencies {
             api(project(":kmath-core"))
-            api("org.jetbrains.kotlinx:kotlinx-coroutines-core-common:${Scientifik.coroutinesVersion}")
-        }
-    }
-    jvmMain {
-        dependencies {
-            api("org.jetbrains.kotlinx:kotlinx-coroutines-core:${Scientifik.coroutinesVersion}")
-        }
-    }
-    jsMain {
-        dependencies {
-            api("org.jetbrains.kotlinx:kotlinx-coroutines-core-js:${Scientifik.coroutinesVersion}")
         }
     }
 }

kmath-functions/src/commonMain/kotlin/scientifik.kmath.functions/Piecewise.kt

@@ -1,2 +0,0 @@
-package scientifik.kmath.functions
-

kmath-functions/src/commonMain/kotlin/scientifik.kmath.functions/Polynomial.kt

@@ -1,4 +1,94 @@
 package scientifik.kmath.functions
 
-interface Polynomial {
+import scientifik.kmath.operations.RealField
+import scientifik.kmath.operations.Ring
+import scientifik.kmath.operations.Space
+import kotlin.jvm.JvmName
+import kotlin.math.max
+import kotlin.math.pow
+
+/**
+ * Polynomial coefficients without fixation on specific context they are applied to
+ * @param coefficients constant is the leftmost coefficient
+ */
+inline class Polynomial<T : Any>(val coefficients: List<T>) {
+    constructor(vararg coefficients: T) : this(coefficients.toList())
 }
+
+fun Polynomial<Double>.value() =
+    coefficients.reduceIndexed { index: Int, acc: Double, d: Double -> acc + d.pow(index) }
+
+
+fun <T : Any, C : Ring<T>> Polynomial<T>.value(ring: C, arg: T): T = ring.run {
+    if( coefficients.isEmpty()) return@run zero
+    var res = coefficients.first()
+    var powerArg = arg
+    for( index in 1 until coefficients.size){
+        res += coefficients[index]*powerArg
+        //recalculating power on each step to avoid power costs on long polynomials
+        powerArg *= arg
+    }
+    return@run res
+}
+
+/**
+ * Represent a polynomial as a context-dependent function
+ */
+fun <T : Any, C : Ring<T>> Polynomial<T>.asMathFunction(): MathFunction<T, out C, T> = object : MathFunction<T, C, T> {
+    override fun C.invoke(arg: T): T = value(this, arg)
+}
+
+/**
+ * Represent the polynomial as a regular context-less function
+ */
+fun <T : Any, C : Ring<T>> Polynomial<T>.asFunction(ring: C): (T) -> T = { value(ring, it) }
+
+@JvmName("asRealUFunction")
+fun Polynomial<Double>.asFunction(): (Double) -> Double = asFunction(RealField)
+
+/**
+ * An algebra for polynomials
+ */
+class PolynomialSpace<T : Any, C : Ring<T>>(val ring: C) : Space<Polynomial<T>> {
+
+    override fun add(a: Polynomial<T>, b: Polynomial<T>): Polynomial<T> {
+        val dim = max(a.coefficients.size, b.coefficients.size)
+        ring.run {
+            return Polynomial(List(dim) { index ->
+                a.coefficients.getOrElse(index) { zero } + b.coefficients.getOrElse(index) { zero }
+            })
+        }
+    }
+
+    override fun multiply(a: Polynomial<T>, k: Number): Polynomial<T> {
+        ring.run {
+            return Polynomial(List(a.coefficients.size) { index -> a.coefficients[index] * k })
+        }
+    }
+
+    override val zero: Polynomial<T> = Polynomial(emptyList())
+
+    operator fun Polynomial<T>.invoke(arg: T): T = value(ring, arg)
+}
+
+fun <T : Any, C : Ring<T>, R> C.polynomial(block: PolynomialSpace<T, C>.() -> R): R {
+    return PolynomialSpace(this).run(block)
+}
+
+class PiecewisePolynomial<T : Comparable<T>> internal constructor(
+    val lowerBoundary: T,
+    val pieces: List<Pair<T, Polynomial<T>>>
+)
+
+private fun <T : Comparable<T>> PiecewisePolynomial<T>.findPiece(arg: T): Polynomial<T>? {
+    if (arg < lowerBoundary || arg > pieces.last().first) return null
+    return pieces.first { arg < it.first }.second
+}
+
+/**
+ * Return a value of polynomial function with given [ring] an given [arg] or null if argument is outside of piecewise definition.
+ */
+fun <T : Comparable<T>, C : Ring<T>> PiecewisePolynomial<T>.value(ring: C, arg: T): T? =
+    findPiece(arg)?.value(ring, arg)
+
+fun <T : Comparable<T>, C : Ring<T>> PiecewisePolynomial<T>.asFunction(ring: C): (T) -> T? = { value(ring, it) }

kmath-functions/src/commonMain/kotlin/scientifik.kmath.functions/functions.kt

@@ -1,54 +1,33 @@
-package scientifik.kmath.misc
+package scientifik.kmath.functions
 
+import scientifik.kmath.operations.Algebra
 import scientifik.kmath.operations.RealField
-import scientifik.kmath.operations.SpaceOperations
-import kotlin.jvm.JvmName
 
 /**
  * A regular function that could be called only inside specific algebra context
+ * @param T source type
+ * @param C source algebra constraint
+ * @param R result type
  */
-interface UFunction<T, C : SpaceOperations<T>>  {
-    operator fun C.invoke(arg: T): T
+interface MathFunction<T, C : Algebra<T>, R> {
+    operator fun C.invoke(arg: T): R
 }
 
+fun <R> MathFunction<Double, RealField, R>.invoke(arg: Double): R = RealField.invoke(arg)
+
 /**
- * A suspendable univariate function defined in algebraic context
+ * A suspendable function defined in algebraic context
  */
-interface USFunction<T, C : SpaceOperations<T>> {
-    suspend operator fun C.invoke(arg: T): T
+interface SuspendableMathFunction<T, C : Algebra<T>, R> {
+    suspend operator fun C.invoke(arg: T): R
 }
 
-suspend fun USFunction<Double, RealField>.invoke(arg: Double) = RealField.invoke(arg)
+suspend fun <R> SuspendableMathFunction<Double, RealField, R>.invoke(arg: Double) = RealField.invoke(arg)
 
 
-interface MFunction<T, C : SpaceOperations<T>> {
-    /**
-     * The input dimension of the function
-     */
-    val dimension: UInt
-
-    operator fun C.invoke(vararg args: T): T
-}
-
 /**
- * A suspendable multivariate (N->1) function defined on algebraic context
+ * A parametric function with parameter
  */
-interface MSFunction<T, C : SpaceOperations<T>> {
-    /**
-     * The input dimension of the function
-     */
-    val dimension: UInt
-
-    suspend operator fun C.invoke(vararg args: T): T
-}
-
-suspend fun MSFunction<Double, RealField>.invoke(args: DoubleArray) = RealField.invoke(*args.toTypedArray())
-@JvmName("varargInvoke")
-suspend fun MSFunction<Double, RealField>.invoke(vararg args: Double) = RealField.invoke(*args.toTypedArray())
-
-/**
- * A suspendable parametric function with parameter
- */
-interface PSFunction<T, P, C : SpaceOperations<T>> {
-    suspend operator fun C.invoke(arg: T, parameter: P): T
+interface ParametricFunction<T, P, C : Algebra<T>> {
+    operator fun C.invoke(arg: T, parameter: P): T
 }

kmath-functions/src/commonMain/kotlin/scientifik/kmath/interpolation/Interpolator.kt

@@ -1,7 +1,21 @@
 package scientifik.kmath.interpolation
 
-import scientifik.kmath.functions.MathFunction
+import scientifik.kmath.functions.PiecewisePolynomial
+import scientifik.kmath.functions.value
+import scientifik.kmath.operations.Ring
 
 interface Interpolator<X, Y> {
-    fun interpolate(points: Collection<Pair<X, Y>>): MathFunction<X, *, Y>
+    fun interpolate(points: Collection<Pair<X, Y>>): (X) -> Y
+}
+
+interface PolynomialInterpolator<T : Comparable<T>> : Interpolator<T, T> {
+    val algebra: Ring<T>
+
+    fun getDefaultValue(): T = error("Out of bounds")
+
+    fun interpolatePolynomials(points: Collection<Pair<T, T>>): PiecewisePolynomial<T>
+
+    override fun interpolate(points: Collection<Pair<T, T>>): (T) -> T = { x ->
+        interpolatePolynomials(points).value(algebra, x) ?: getDefaultValue()
+    }
 }

kmath-functions/src/commonMain/kotlin/scientifik/kmath/interpolation/LinearInterpolator.kt

@@ -1,4 +1,24 @@
 package scientifik.kmath.interpolation
 
-class LinearInterpolator {
+import scientifik.kmath.functions.PiecewisePolynomial
+import scientifik.kmath.functions.Polynomial
+import scientifik.kmath.operations.Field
+
+/**
+ * Reference JVM implementation: https://github.com/apache/commons-math/blob/master/src/main/java/org/apache/commons/math4/analysis/interpolation/LinearInterpolator.java
+ */
+class LinearInterpolator<T : Comparable<T>>(override val algebra: Field<T>) : PolynomialInterpolator<T> {
+
+    override fun interpolatePolynomials(points: Collection<Pair<T, T>>): PiecewisePolynomial<T> = algebra.run {
+        //sorting points
+        val sorted = points.sortedBy { it.first }
+
+        val pairs: List<Pair<T, Polynomial<T>>> = (0 until points.size - 1).map { i ->
+            val slope = (sorted[i + 1].second - sorted[i].second) / (sorted[i + 1].first - sorted[i].first)
+            val const = sorted[i].second - slope * sorted[i].first
+            sorted[i + 1].first to Polynomial(const, slope)
+        }
+
+        return PiecewisePolynomial(sorted.first().first, pairs)
+    }
 }

kmath-functions/src/commonTest/kotlin/scientifik/kmath/interpolation/LinearInterpolatorTest.kt

@@ -1,5 +1,26 @@
 package scientifik.kmath.interpolation
 
-import org.junit.Assert.*
+import scientifik.kmath.functions.asFunction
+import scientifik.kmath.operations.RealField
+import kotlin.test.Test
+import kotlin.test.assertEquals
 
-class LinearInterpolatorTest
+
+class LinearInterpolatorTest {
+    @Test
+    fun testInterpolation() {
+        val data = listOf(
+            0.0 to 0.0,
+            1.0 to 1.0,
+            2.0 to 3.0,
+            3.0 to 4.0
+        )
+        val polynomial = LinearInterpolator(RealField).interpolatePolynomials(data)
+        val function = polynomial.asFunction(RealField)
+
+//        assertEquals(null, function(-1.0))
+//        assertEquals(0.5, function(0.5))
+        assertEquals(2.0, function(1.5))
+        assertEquals(3.0, function(2.0))
+    }
+}

settings.gradle.kts

@@ -1,10 +1,10 @@
 pluginManagement {
 
     plugins {
-        id("scientifik.mpp") version "0.2.5"
-        id("scientifik.jvm") version "0.2.5"
-        id("scientifik.atomic") version "0.2.5"
-        id("scientifik.publish") version "0.2.5"
+        id("scientifik.mpp") version "0.3.1"
+        id("scientifik.jvm") version "0.3.1"
+        id("scientifik.atomic") version "0.3.1"
+        id("scientifik.publish") version "0.3.1"
     }
 
     repositories {
@@ -19,7 +19,7 @@ pluginManagement {
     resolutionStrategy {
         eachPlugin {
             when (requested.id.id) {
-                "scientifik.mpp", "scientifik.publish" -> useModule("scientifik:gradle-tools:${requested.version}")
+                "scientifik.mpp", "scientifik.jvm", "scientifik.publish" -> useModule("scientifik:gradle-tools:${requested.version}")
             }
         }
     }
@@ -29,6 +29,7 @@ rootProject.name = "kmath"
 include(
     ":kmath-memory",
     ":kmath-core",
+    ":kmath-functions",
 //    ":kmath-io",
     ":kmath-coroutines",
     ":kmath-histograms",