Feature: Polynomials and rational functions #469

Merged
lounres merged 132 commits from feature/polynomials into dev 2022-07-28 18:04:06 +03:00
3 changed files with 58 additions and 49 deletions
Showing only changes of commit 0a5122a974 - Show all commits

View File

@ -5,6 +5,7 @@
package space.kscience.kmath.functions
import space.kscience.kmath.operations.invoke
import space.kscience.kmath.operations.Ring
import space.kscience.kmath.operations.ScaleOperations
import kotlin.contracts.InvocationKind
@ -132,7 +133,59 @@ public fun <C, A: Ring<C>> NumberedPolynomialSpace<C, A>.NumberedPolynomial(vara
//context(NumberedPolynomialSpace<C, A>)
//public fun <C, A: Ring<C>> Symbol.asNumberedPolynomial() : NumberedPolynomial<C> = NumberedPolynomial<C>(mapOf(mapOf(this to 1u) to constantOne))
public fun <C> C.asNumberedPolynomial() : NumberedPolynomial<C> = NumberedPolynomial<C>(mapOf(emptyList<UInt>() to this))
//context(A)
//public fun <C> C.asNumberedPolynomial() : NumberedPolynomial<C> = NumberedPolynomial<C>(mapOf(emptyList<UInt>() to this))
@DslMarker
internal annotation class NumberedPolynomialConstructorDSL
@NumberedPolynomialConstructorDSL
public class NumberedPolynomialTermSignatureBuilder {
private val signature: MutableList<UInt> = ArrayList()
public fun build(): List<UInt> = signature
public infix fun Int.inPowerOf(deg: UInt) {
if (this > signature.lastIndex) {
signature.addAll(List(this - signature.lastIndex - 1) { 0u })
signature.add(deg)
} else {
signature[this] = deg
}
}
public infix fun Int.to(deg: UInt): Unit = this inPowerOf deg
}
@NumberedPolynomialConstructorDSL
public class NumberedPolynomialBuilderOverRing<C> internal constructor(internal val context: Ring<C>, capacity: Int = 0) {
private val coefficients: MutableMap<List<UInt>, C> = LinkedHashMap(capacity)
public fun build(): NumberedPolynomial<C> = NumberedPolynomial<C>(coefficients)
public operator fun C.invoke(block: NumberedPolynomialTermSignatureBuilder.() -> Unit) {
val signature = NumberedPolynomialTermSignatureBuilder().apply(block).build()
coefficients[signature] = context { coefficients.getOrElse(signature) { zero } + this@invoke }
}
}
@NumberedPolynomialConstructorDSL
public class NumberedPolynomialBuilderOverPolynomialSpace<C> internal constructor(internal val context: NumberedPolynomialSpace<C, *>, capacity: Int = 0) {
private val coefficients: MutableMap<List<UInt>, C> = LinkedHashMap(capacity)
public fun build(): NumberedPolynomial<C> = NumberedPolynomial<C>(coefficients)
public operator fun C.invoke(block: NumberedPolynomialTermSignatureBuilder.() -> Unit) {
val signature = NumberedPolynomialTermSignatureBuilder().apply(block).build()
coefficients[signature] = context { coefficients.getOrElse(signature) { constantZero } + this@invoke }
}
}
@NumberedPolynomialConstructorDSL
@Suppress("FunctionName")
public fun <C, A: Ring<C>> A.NumberedPolynomial(block: NumberedPolynomialBuilderOverRing<C>.() -> Unit) : NumberedPolynomial<C> = NumberedPolynomialBuilderOverRing(this).apply(block).build()
@NumberedPolynomialConstructorDSL
@Suppress("FunctionName")
public fun <C, A: Ring<C>> A.NumberedPolynomial(capacity: Int, block: NumberedPolynomialBuilderOverRing<C>.() -> Unit) : NumberedPolynomial<C> = NumberedPolynomialBuilderOverRing(this, capacity).apply(block).build()
@NumberedPolynomialConstructorDSL
@Suppress("FunctionName")
public fun <C, A: Ring<C>> NumberedPolynomialSpace<C, A>.NumberedPolynomial(block: NumberedPolynomialBuilderOverPolynomialSpace<C>.() -> Unit) : NumberedPolynomial<C> = NumberedPolynomialBuilderOverPolynomialSpace(this).apply(block).build()
@NumberedPolynomialConstructorDSL
@Suppress("FunctionName")
public fun <C, A: Ring<C>> NumberedPolynomialSpace<C, A>.NumberedPolynomial(capacity: Int, block: NumberedPolynomialBuilderOverPolynomialSpace<C>.() -> Unit) : NumberedPolynomial<C> = NumberedPolynomialBuilderOverPolynomialSpace(this, capacity).apply(block).build()
/**
* Space of polynomials.
@ -604,4 +657,7 @@ public open class NumberedPolynomialSpace<C, A : Ring<C>>(
for ((degs, c) in this) if (c.isZero()) this.remove(degs)
}
}
// TODO: Move to other constructors with context receiver
public fun C.asNumberedPolynomial() : NumberedPolynomial<C> = NumberedPolynomial<C>(mapOf(emptyList<UInt>() to this))
}

View File

@ -17,53 +17,6 @@ import kotlin.contracts.contract
// TODO: Docs
//// TODO: Reuse underlying ring extensions
//
//context(LabeledPolynomialSpace<C, A>)
//@Suppress("NOTHING_TO_INLINE")
//fun <C, A: Ring<C>> numberConstant(value: Int): C = ring { number<C>(value) }
//
//context(LabeledPolynomialSpace<C, A>)
//fun <C, A: Ring<C>> power(arg: C, pow: UInt): C = ring { power(arg, pow) }
//
//context(LabeledPolynomialSpace<C, A>)
//fun <C, A: Ring<C>> multiplyWithPower(base: C, arg: C, pow: UInt): C = ring { multiplyWithPower<C>(base, arg, pow) }
//
//
//context(LabeledPolynomialSpace<C, A>)
//fun <C, A: Ring<C>> power(arg: Symbol, pow: UInt): LabeledPolynomial<C> =
// if (pow == 0U) one
// else LabeledPolynomial<C>(mapOf(
// mapOf(arg to pow) to constantOne
// ))
//
//
//context(LabeledPolynomialSpace<C, A>)
//fun <C, A: Ring<C>> number(value: Int): LabeledPolynomial<C> = ring { LabeledPolynomial<C>(mapOf(emptyMap<Symbol, UInt>() to number<C>(value))) }
//
//context(LabeledPolynomialSpace<C, A>)
//fun <C, A: Ring<C>> multiplyWithPower(base: LabeledPolynomial<C>, arg: LabeledPolynomial<C>, pow: UInt): LabeledPolynomial<C> =
// when {
// arg.isZero() && pow > 0U -> base
// arg.isOne() -> base
// arg.isMinusOne() -> if (pow % 2U == 0U) base else -base
// else -> multiplyWithPowerInternalLogic(base, arg, pow)
// }
//
//// Trivial but very slow
//context(LabeledPolynomialSpace<C, A>)
//internal tailrec fun <C, A: Ring<C>> multiplyWithPowerInternalLogic(base: LabeledPolynomial<C>, arg: LabeledPolynomial<C>, exponent: UInt): LabeledPolynomial<C> =
// when {
// exponent == 0U -> base
// exponent == 1U -> base * arg
// exponent % 2U == 0U -> multiplyWithPowerInternalLogic(base, arg * arg, exponent / 2U)
// exponent % 2U == 1U -> multiplyWithPowerInternalLogic(base * arg, arg * arg, exponent / 2U)
// else -> error("Error in raising ring instant by unsigned integer: got reminder by division by 2 different from 0 and 1")
// }
//
/**
* Creates a [LabeledPolynomialSpace] over a received ring.
*/

View File

@ -124,7 +124,7 @@ public fun <C> ListPolynomial<C>.substitute(ring: Ring<C>, arg: ListPolynomial<C
val thisDegree = coefficients.indexOfLast { it != zero }
if (thisDegree == -1) return ListPolynomial(emptyList())
val argDegree = arg.coefficients.indexOfLast { it != zero }
if (argDegree == -1) return coefficients[0].asPolynomial()
if (argDegree == -1) return coefficients[0].asListPolynomial()
val constantZero = zero
val resultCoefs: MutableList<C> = MutableList(thisDegree * argDegree + 1) { constantZero }
resultCoefs[0] = coefficients[thisDegree]