kscience/kmath
9
0
Fork 1
Code Issues Pull Requests 6 Packages Projects Releases 16 Wiki Activity

Feature: Polynomials and rational functions #469

Merged
lounres merged 132 commits from feature/polynomials into dev 2022-07-28 18:04:06 +03:00
Conversation 19 Commits 132 Files Changed 51 +20 -358
2 changed files with 20 additions and 358 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files
Showing only changes of commit bb5e638b31 - Show all commits

189
kmath-functions/src/commonMain/kotlin/space/kscience/kmath/functions/labeledPolynomialUtil.kt
View File

@ -54,7 +54,7 @@ import kotlin.contracts.*
// else -> multiplyWithPowerInternalLogic(base, arg, pow)
// }
//
//// Trivial but slow as duck
//// 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 {
@ -71,7 +71,15 @@ import kotlin.contracts.*
// region Utilities
// TODO: Docs
/**
* Crates a [LabeledPolynomialSpace] over received ring.
*/
public fun <C, A : Ring<C>> A.labeledPolynomial(): LabeledPolynomialSpace<C, A> =
LabeledPolynomialSpace(this)
/**
* Crates a [LabeledPolynomialSpace]'s scope over received ring.
*/
@OptIn(ExperimentalContracts::class)
public inline fun <C, A : Ring<C>, R> A.labeledPolynomial(block: LabeledPolynomialSpace<C, A>.() -> R): R {
contract { callsInPlace(block, InvocationKind.EXACTLY_ONCE) }
@ -254,64 +262,6 @@ public inline fun <C, A : Ring<C>, R> A.labeledPolynomial(block: LabeledPolynomi
//// region Field case
//
//operator fun <T: Field<T>> Polynomial<T>.div(other: Polynomial<T>): Polynomial<T> {
// if (other.isZero()) throw ArithmeticException("/ by zero")
// if (isZero()) return this
//
// fun Map<List<Int>, T>.leadingTerm() =
// this
// .asSequence()
// .map { Pair(it.key, it.value) }
// .reduce { (accDegs, accC), (listDegs, listC) ->
// for (i in 0..accDegs.lastIndex) {
// if (accDegs[i] > listDegs.getOrElse(i) { 0 }) return@reduce accDegs to accC
// if (accDegs[i] < listDegs.getOrElse(i) { 0 }) return@reduce listDegs to listC
// }
// if (accDegs.size < listDegs.size) listDegs to listC else accDegs to accC
// }
//
// var thisCoefficients = coefficients.toMutableMap()
// val otherCoefficients = other.coefficients
// val quotientCoefficients = HashMap<List<Int>, T>()
//
// var (thisLeadingTermDegs, thisLeadingTermC) = thisCoefficients.leadingTerm()
// val (otherLeadingTermDegs, otherLeadingTermC) = otherCoefficients.leadingTerm()
//
// while (
// thisLeadingTermDegs.size >= otherLeadingTermDegs.size &&
// (0..otherLeadingTermDegs.lastIndex).all { thisLeadingTermDegs[it] >= otherLeadingTermDegs[it] }
// ) {
// val multiplierDegs =
// thisLeadingTermDegs
// .mapIndexed { index, deg -> deg - otherLeadingTermDegs.getOrElse(index) { 0 } }
// .cleanUp()
// val multiplierC = thisLeadingTermC / otherLeadingTermC
//
// quotientCoefficients[multiplierDegs] = multiplierC
//
// for ((degs, t) in otherCoefficients) {
// val productDegs =
// (0..max(degs.lastIndex, multiplierDegs.lastIndex))
// .map { degs.getOrElse(it) { 0 } + multiplierDegs.getOrElse(it) { 0 } }
// .cleanUp()
// val productC = t * multiplierC
// thisCoefficients[productDegs] =
// if (productDegs in thisCoefficients) thisCoefficients[productDegs]!! - productC else -productC
// }
//
// thisCoefficients = thisCoefficients.filterValues { it.isNotZero() }.toMutableMap()
//
// if (thisCoefficients.isEmpty())
// return Polynomial(quotientCoefficients, toCheckInput = false)
//
// val t = thisCoefficients.leadingTerm()
// thisLeadingTermDegs = t.first
// thisLeadingTermC = t.second
// }
//
// return Polynomial(quotientCoefficients, toCheckInput = false)
//}
//
//operator fun <T: Field<T>> Polynomial<T>.div(other: T): Polynomial<T> =
// if (other.isZero()) throw ArithmeticException("/ by zero")
// else
@ -321,125 +271,6 @@ public inline fun <C, A : Ring<C>, R> A.labeledPolynomial(block: LabeledPolynomi
// toCheckInput = false
// )
//
//operator fun <T: Field<T>> Polynomial<T>.rem(other: Polynomial<T>): Polynomial<T> {
// if (other.isZero()) throw ArithmeticException("/ by zero")
// if (isZero()) return this
//
// fun Map<List<Int>, T>.leadingTerm() =
// this
// .asSequence()
// .map { Pair(it.key, it.value) }
// .reduce { (accDegs, accC), (listDegs, listC) ->
// for (i in 0..accDegs.lastIndex) {
// if (accDegs[i] > listDegs.getOrElse(i) { 0 }) return@reduce accDegs to accC
// if (accDegs[i] < listDegs.getOrElse(i) { 0 }) return@reduce listDegs to listC
// }
// if (accDegs.size < listDegs.size) listDegs to listC else accDegs to accC
// }
//
// var thisCoefficients = coefficients.toMutableMap()
// val otherCoefficients = other.coefficients
//
// var (thisLeadingTermDegs, thisLeadingTermC) = thisCoefficients.leadingTerm()
// val (otherLeadingTermDegs, otherLeadingTermC) = otherCoefficients.leadingTerm()
//
// while (
// thisLeadingTermDegs.size >= otherLeadingTermDegs.size &&
// (0..otherLeadingTermDegs.lastIndex).all { thisLeadingTermDegs[it] >= otherLeadingTermDegs[it] }
// ) {
// val multiplierDegs =
// thisLeadingTermDegs
// .mapIndexed { index, deg -> deg - otherLeadingTermDegs.getOrElse(index) { 0 } }
// .cleanUp()
// val multiplierC = thisLeadingTermC / otherLeadingTermC
//
// for ((degs, t) in otherCoefficients) {
// val productDegs =
// (0..max(degs.lastIndex, multiplierDegs.lastIndex))
// .map { degs.getOrElse(it) { 0 } + multiplierDegs.getOrElse(it) { 0 } }
// .cleanUp()
// val productC = t * multiplierC
// thisCoefficients[productDegs] =
// if (productDegs in thisCoefficients) thisCoefficients[productDegs]!! - productC else -productC
// }
//
// thisCoefficients = thisCoefficients.filterValues { it.isNotZero() }.toMutableMap()
//
// if (thisCoefficients.isEmpty())
// return Polynomial(thisCoefficients, toCheckInput = false)
//
// val t = thisCoefficients.leadingTerm()
// thisLeadingTermDegs = t.first
// thisLeadingTermC = t.second
// }
//
// return Polynomial(thisCoefficients, toCheckInput = false)
//}
//
//infix fun <T: Field<T>> Polynomial<T>.divrem(other: Polynomial<T>): Polynomial.Companion.DividingResult<T> {
// if (other.isZero()) throw ArithmeticException("/ by zero")
// if (isZero()) return Polynomial.Companion.DividingResult(this, this)
//
// fun Map<List<Int>, T>.leadingTerm() =
// this
// .asSequence()
// .map { Pair(it.key, it.value) }
// .reduce { (accDegs, accC), (listDegs, listC) ->
// for (i in 0..accDegs.lastIndex) {
// if (accDegs[i] > listDegs.getOrElse(i) { 0 }) return@reduce accDegs to accC
// if (accDegs[i] < listDegs.getOrElse(i) { 0 }) return@reduce listDegs to listC
// }
// if (accDegs.size < listDegs.size) listDegs to listC else accDegs to accC
// }
//
// var thisCoefficients = coefficients.toMutableMap()
// val otherCoefficients = other.coefficients
// val quotientCoefficients = HashMap<List<Int>, T>()
//
// var (thisLeadingTermDegs, thisLeadingTermC) = thisCoefficients.leadingTerm()
// val (otherLeadingTermDegs, otherLeadingTermC) = otherCoefficients.leadingTerm()
//
// while (
// thisLeadingTermDegs.size >= otherLeadingTermDegs.size &&
// (0..otherLeadingTermDegs.lastIndex).all { thisLeadingTermDegs[it] >= otherLeadingTermDegs[it] }
// ) {
// val multiplierDegs =
// thisLeadingTermDegs
// .mapIndexed { index, deg -> deg - otherLeadingTermDegs.getOrElse(index) { 0 } }
// .cleanUp()
// val multiplierC = thisLeadingTermC / otherLeadingTermC
//
// quotientCoefficients[multiplierDegs] = multiplierC
//
// for ((degs, t) in otherCoefficients) {
// val productDegs =
// (0..max(degs.lastIndex, multiplierDegs.lastIndex))
// .map { degs.getOrElse(it) { 0 } + multiplierDegs.getOrElse(it) { 0 } }
// .cleanUp()
// val productC = t * multiplierC
// thisCoefficients[productDegs] =
// if (productDegs in thisCoefficients) thisCoefficients[productDegs]!! - productC else -productC
// }
//
// thisCoefficients = thisCoefficients.filterValues { it.isNotZero() }.toMutableMap()
//
// if (thisCoefficients.isEmpty())
// return Polynomial.Companion.DividingResult(
// Polynomial(quotientCoefficients, toCheckInput = false),
// Polynomial(thisCoefficients, toCheckInput = false)
// )
//
// val t = thisCoefficients.leadingTerm()
// thisLeadingTermDegs = t.first
// thisLeadingTermC = t.second
// }
//
// return Polynomial.Companion.DividingResult(
// Polynomial(quotientCoefficients, toCheckInput = false),
// Polynomial(thisCoefficients, toCheckInput = false)
// )
//}
//
//// endregion
// endregion

189
kmath-functions/src/commonMain/kotlin/space/kscience/kmath/functions/numberedPolynomialUtil.kt
View File

@ -38,7 +38,7 @@ import kotlin.jvm.JvmName
// else -> multiplyWithPowerInternalLogic(base, arg, pow)
// }
//
//// Trivial but slow as duck
//// Trivial but very slow
//context(NumberedPolynomialSpace<C, A>)
//internal tailrec fun <C, A: Ring<C>> multiplyWithPowerInternalLogic(base: NumberedPolynomial<C>, arg: NumberedPolynomial<C>, exponent: UInt): NumberedPolynomial<C> =
// when {
@ -55,7 +55,15 @@ import kotlin.jvm.JvmName
// region Utilities
// TODO: Docs
/**
* Crates a [NumberedPolynomialSpace] over received ring.
*/
public fun <C, A : Ring<C>> A.numberedPolynomial(): NumberedPolynomialSpace<C, A> =
NumberedPolynomialSpace(this)
/**
* Crates a [NumberedPolynomialSpace]'s scope over received ring.
*/
@OptIn(ExperimentalContracts::class)
public inline fun <C, A : Ring<C>, R> A.numberedPolynomial(block: NumberedPolynomialSpace<C, A>.() -> R): R {
contract { callsInPlace(block, InvocationKind.EXACTLY_ONCE) }
@ -282,64 +290,6 @@ public inline fun <C, A : Ring<C>, R> A.numberedPolynomial(block: NumberedPolyno
//// region Field case
//
//operator fun <T: Field<T>> Polynomial<T>.div(other: Polynomial<T>): Polynomial<T> {
// if (other.isZero()) throw ArithmeticException("/ by zero")
// if (isZero()) return this
//
// fun Map<List<Int>, T>.leadingTerm() =
// this
// .asSequence()
// .map { Pair(it.key, it.value) }
// .reduce { (accDegs, accC), (listDegs, listC) ->
// for (i in 0..accDegs.lastIndex) {
// if (accDegs[i] > listDegs.getOrElse(i) { 0 }) return@reduce accDegs to accC
// if (accDegs[i] < listDegs.getOrElse(i) { 0 }) return@reduce listDegs to listC
// }
// if (accDegs.size < listDegs.size) listDegs to listC else accDegs to accC
// }
//
// var thisCoefficients = coefficients.toMutableMap()
// val otherCoefficients = other.coefficients
// val quotientCoefficients = HashMap<List<Int>, T>()
//
// var (thisLeadingTermDegs, thisLeadingTermC) = thisCoefficients.leadingTerm()
// val (otherLeadingTermDegs, otherLeadingTermC) = otherCoefficients.leadingTerm()
//
// while (
// thisLeadingTermDegs.size >= otherLeadingTermDegs.size &&
// (0..otherLeadingTermDegs.lastIndex).all { thisLeadingTermDegs[it] >= otherLeadingTermDegs[it] }
// ) {
// val multiplierDegs =
// thisLeadingTermDegs
// .mapIndexed { index, deg -> deg - otherLeadingTermDegs.getOrElse(index) { 0 } }
// .cleanUp()
// val multiplierC = thisLeadingTermC / otherLeadingTermC
//
// quotientCoefficients[multiplierDegs] = multiplierC
//
// for ((degs, t) in otherCoefficients) {
// val productDegs =
// (0..max(degs.lastIndex, multiplierDegs.lastIndex))
// .map { degs.getOrElse(it) { 0 } + multiplierDegs.getOrElse(it) { 0 } }
// .cleanUp()
// val productC = t * multiplierC
// thisCoefficients[productDegs] =
// if (productDegs in thisCoefficients) thisCoefficients[productDegs]!! - productC else -productC
// }
//
// thisCoefficients = thisCoefficients.filterValues { it.isNotZero() }.toMutableMap()
//
// if (thisCoefficients.isEmpty())
// return Polynomial(quotientCoefficients, toCheckInput = false)
//
// val t = thisCoefficients.leadingTerm()
// thisLeadingTermDegs = t.first
// thisLeadingTermC = t.second
// }
//
// return Polynomial(quotientCoefficients, toCheckInput = false)
//}
//
//operator fun <T: Field<T>> Polynomial<T>.div(other: T): Polynomial<T> =
// if (other.isZero()) throw ArithmeticException("/ by zero")
// else
@ -349,125 +299,6 @@ public inline fun <C, A : Ring<C>, R> A.numberedPolynomial(block: NumberedPolyno
// toCheckInput = false
// )
//
//operator fun <T: Field<T>> Polynomial<T>.rem(other: Polynomial<T>): Polynomial<T> {
// if (other.isZero()) throw ArithmeticException("/ by zero")
// if (isZero()) return this
//
// fun Map<List<Int>, T>.leadingTerm() =
// this
// .asSequence()
// .map { Pair(it.key, it.value) }
// .reduce { (accDegs, accC), (listDegs, listC) ->
// for (i in 0..accDegs.lastIndex) {
// if (accDegs[i] > listDegs.getOrElse(i) { 0 }) return@reduce accDegs to accC
// if (accDegs[i] < listDegs.getOrElse(i) { 0 }) return@reduce listDegs to listC
// }
// if (accDegs.size < listDegs.size) listDegs to listC else accDegs to accC
// }
//
// var thisCoefficients = coefficients.toMutableMap()
// val otherCoefficients = other.coefficients
//
// var (thisLeadingTermDegs, thisLeadingTermC) = thisCoefficients.leadingTerm()
// val (otherLeadingTermDegs, otherLeadingTermC) = otherCoefficients.leadingTerm()
//
// while (
// thisLeadingTermDegs.size >= otherLeadingTermDegs.size &&
// (0..otherLeadingTermDegs.lastIndex).all { thisLeadingTermDegs[it] >= otherLeadingTermDegs[it] }
// ) {
// val multiplierDegs =
// thisLeadingTermDegs
// .mapIndexed { index, deg -> deg - otherLeadingTermDegs.getOrElse(index) { 0 } }
// .cleanUp()
// val multiplierC = thisLeadingTermC / otherLeadingTermC
//
// for ((degs, t) in otherCoefficients) {
// val productDegs =
// (0..max(degs.lastIndex, multiplierDegs.lastIndex))
// .map { degs.getOrElse(it) { 0 } + multiplierDegs.getOrElse(it) { 0 } }
// .cleanUp()
// val productC = t * multiplierC
// thisCoefficients[productDegs] =
// if (productDegs in thisCoefficients) thisCoefficients[productDegs]!! - productC else -productC
// }
//
// thisCoefficients = thisCoefficients.filterValues { it.isNotZero() }.toMutableMap()
//
// if (thisCoefficients.isEmpty())
// return Polynomial(thisCoefficients, toCheckInput = false)
//
// val t = thisCoefficients.leadingTerm()
// thisLeadingTermDegs = t.first
// thisLeadingTermC = t.second
// }
//
// return Polynomial(thisCoefficients, toCheckInput = false)
//}
//
//infix fun <T: Field<T>> Polynomial<T>.divrem(other: Polynomial<T>): Polynomial.Companion.DividingResult<T> {
// if (other.isZero()) throw ArithmeticException("/ by zero")
// if (isZero()) return Polynomial.Companion.DividingResult(this, this)
//
// fun Map<List<Int>, T>.leadingTerm() =
// this
// .asSequence()
// .map { Pair(it.key, it.value) }
// .reduce { (accDegs, accC), (listDegs, listC) ->
// for (i in 0..accDegs.lastIndex) {
// if (accDegs[i] > listDegs.getOrElse(i) { 0 }) return@reduce accDegs to accC
// if (accDegs[i] < listDegs.getOrElse(i) { 0 }) return@reduce listDegs to listC
// }
// if (accDegs.size < listDegs.size) listDegs to listC else accDegs to accC
// }
//
// var thisCoefficients = coefficients.toMutableMap()
// val otherCoefficients = other.coefficients
// val quotientCoefficients = HashMap<List<Int>, T>()
//
// var (thisLeadingTermDegs, thisLeadingTermC) = thisCoefficients.leadingTerm()
// val (otherLeadingTermDegs, otherLeadingTermC) = otherCoefficients.leadingTerm()
//
// while (
// thisLeadingTermDegs.size >= otherLeadingTermDegs.size &&
// (0..otherLeadingTermDegs.lastIndex).all { thisLeadingTermDegs[it] >= otherLeadingTermDegs[it] }
// ) {
// val multiplierDegs =
// thisLeadingTermDegs
// .mapIndexed { index, deg -> deg - otherLeadingTermDegs.getOrElse(index) { 0 } }
// .cleanUp()
// val multiplierC = thisLeadingTermC / otherLeadingTermC
//
// quotientCoefficients[multiplierDegs] = multiplierC
//
// for ((degs, t) in otherCoefficients) {
// val productDegs =
// (0..max(degs.lastIndex, multiplierDegs.lastIndex))
// .map { degs.getOrElse(it) { 0 } + multiplierDegs.getOrElse(it) { 0 } }
// .cleanUp()
// val productC = t * multiplierC
// thisCoefficients[productDegs] =
// if (productDegs in thisCoefficients) thisCoefficients[productDegs]!! - productC else -productC
// }
//
// thisCoefficients = thisCoefficients.filterValues { it.isNotZero() }.toMutableMap()
//
// if (thisCoefficients.isEmpty())
// return Polynomial.Companion.DividingResult(
// Polynomial(quotientCoefficients, toCheckInput = false),
// Polynomial(thisCoefficients, toCheckInput = false)
// )
//
// val t = thisCoefficients.leadingTerm()
// thisLeadingTermDegs = t.first
// thisLeadingTermC = t.second
// }
//
// return Polynomial.Companion.DividingResult(
// Polynomial(quotientCoefficients, toCheckInput = false),
// Polynomial(thisCoefficients, toCheckInput = false)
// )
//}
//
//// endregion
// endregion