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

Added derivative-like functions to Polynomial

Browse Source
This commit is contained in:
Gleb Minaev 2022-03-15 18:10:11 +03:00
parent 1754ae0695
commit 91c9ea61da
2 changed files with 39 additions and 25 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

16
kmath-functions/src/commonMain/kotlin/space/kscience/kmath/functions/Polynomial.kt
View File

@ -14,7 +14,7 @@ import kotlin.math.min
*
* @param coefficients constant is the leftmost coefficient.
*/
public data class Polynomial<T>(public val coefficients: List<T>) : AbstractPolynomial<T> {
public data class Polynomial<C>(public val coefficients: List<C>) : AbstractPolynomial<C> {
override fun toString(): String = "Polynomial$coefficients"
}
@ -44,7 +44,7 @@ internal fun polynomialError(message: Any): Nothing = throw PolynomialError(mess
* [reverse] parameter is true.
*/
@Suppress("FunctionName")
public fun <T> Polynomial(coefficients: List<T>, reverse: Boolean = false): Polynomial<T> =
public fun <C> Polynomial(coefficients: List<C>, reverse: Boolean = false): Polynomial<C> =
Polynomial(with(coefficients) { if (reverse) reversed() else this })
/**
@ -52,10 +52,10 @@ public fun <T> Polynomial(coefficients: List<T>, reverse: Boolean = false): Poly
* [reverse] parameter is true.
*/
@Suppress("FunctionName")
public fun <T> Polynomial(vararg coefficients: T, reverse: Boolean = false): Polynomial<T> =
public fun <C> Polynomial(vararg coefficients: C, reverse: Boolean = false): Polynomial<C> =
Polynomial(with(coefficients) { if (reverse) reversed() else toList() })
public fun <T> T.asPolynomial() : Polynomial<T> = Polynomial(listOf(this))
public fun <C> C.asPolynomial() : Polynomial<C> = Polynomial(listOf(this))
// endregion
@ -352,7 +352,7 @@ public open class PolynomialSpace<C, A : Ring<C>>(
Polynomial(
(0..(thisDegree + otherDegree))
.map { d ->
(max(0, d - otherDegree)..(min(thisDegree, d)))
(max(0, d - otherDegree)..min(thisDegree, d))
.map { coefficients[it] * other.coefficients[d - it] }
.reduce { acc, rational -> acc + rational }
}
@ -370,16 +370,14 @@ public open class PolynomialSpace<C, A : Ring<C>>(
*/
public override fun Polynomial<C>.isOne(): Boolean =
with(coefficients) {
isNotEmpty() &&
asSequence().withIndex().any { (index, c) -> if (index == 0) c.isOne() else c.isZero() } // TODO: It's better to write new methods like `anyIndexed`. But what's better way to do it?
isNotEmpty() && withIndex().any { (index, c) -> if (index == 0) c.isOne() else c.isZero() }
}
/**
* Check if the instant is minus unit polynomial.
*/
public override fun Polynomial<C>.isMinusOne(): Boolean =
with(coefficients) {
isNotEmpty() &&
asSequence().withIndex().any { (index, c) -> if (index == 0) c.isMinusOne() else c.isZero() } // TODO: It's better to write new methods like `anyIndexed`. But what's better way to do it?
isNotEmpty() && withIndex().any { (index, c) -> if (index == 0) c.isMinusOne() else c.isZero() }
}
/**

48
kmath-functions/src/commonMain/kotlin/space/kscience/kmath/functions/polynomialUtil.kt
View File

@ -9,7 +9,8 @@ import space.kscience.kmath.misc.UnstableKMathAPI
import space.kscience.kmath.operations.*
import kotlin.contracts.InvocationKind
import kotlin.contracts.contract
import kotlin.jvm.JvmName
import kotlin.math.max
import kotlin.math.min
import kotlin.math.pow
@ -77,15 +78,32 @@ public fun <C> Polynomial<C>.substitute(ring: Ring<C>, arg: C): C = ring {
return result
}
// TODO: (Waiting for hero) Replace with optimisation: the [result] may be unboxed, and all operations may be performed
// as soon as possible on it
public fun <C> Polynomial<C>.substitute(ring: Ring<C>, arg: Polynomial<C>) : Polynomial<C> = ring.polynomial {
if (coefficients.isEmpty()) return zero
var result: Polynomial<C> = coefficients.last().asPolynomial()
for (j in coefficients.size - 2 downTo 0) {
result = (arg * result) + coefficients[j]
val thisDegree = degree
if (thisDegree == -1) return zero
val argDegree = arg.degree
if (argDegree == -1) return coefficients[0].asPolynomial()
val constantZero = constantZero
val resultCoefs: MutableList<C> = MutableList(thisDegree + argDegree + 1) { constantZero }
val resultCoefsUpdate: MutableList<C> = MutableList(thisDegree + argDegree + 1) { constantZero }
var resultDegree = 0
for (deg in thisDegree downTo 0) {
resultCoefsUpdate[0] = coefficients[deg]
for (updateDeg in 0 .. resultDegree + argDegree) {
var newC = resultCoefsUpdate[updateDeg]
for (deg1 in max(0, updateDeg - argDegree)..min(resultDegree, updateDeg))
newC += resultCoefs[deg1] * arg.coefficients[updateDeg - deg1]
resultCoefsUpdate[updateDeg] = newC
}
resultDegree += argDegree
for (updateDeg in 0 .. resultDegree + argDegree) {
resultCoefs[updateDeg] = resultCoefsUpdate[updateDeg]
resultCoefsUpdate[updateDeg] = constantZero
}
}
return result
return Polynomial<C>(resultCoefs)
}
/**
@ -109,7 +127,7 @@ public fun <C, A : Ring<C>> Polynomial<C>.asPolynomialFunctionOver(ring: A): (Po
public fun <C, A> Polynomial<C>.derivative(
algebra: A,
): Polynomial<C> where A : Ring<C>, A : NumericAlgebra<C> = algebra {
Polynomial(coefficients.drop(1).mapIndexed { index, c -> number(index) * c })
Polynomial(coefficients.drop(1).mapIndexed { index, c -> number(index + 1) * c })
}
/**
@ -120,8 +138,7 @@ public fun <C, A> Polynomial<C>.nthDerivative(
algebra: A,
order: UInt,
): Polynomial<C> where A : Ring<C>, A : NumericAlgebra<C> = algebra {
TODO()
Polynomial(coefficients.drop(order.toInt()).mapIndexed { index, c -> number(index) * c })
Polynomial(coefficients.drop(order.toInt()).mapIndexed { index, c -> (1..order.toInt()).fold(c) { acc, i -> acc * number(index + i) } })
}
/**
@ -133,7 +150,7 @@ public fun <C, A> Polynomial<C>.antiderivative(
): Polynomial<C> where A : Field<C>, A : NumericAlgebra<C> = algebra {
val integratedCoefficients = buildList(coefficients.size + 1) {
add(zero)
coefficients.forEachIndexed{ index, t -> add(t / number(index + 1)) }
coefficients.mapIndexedTo(this) { index, t -> t / number(index + 1) }
}
Polynomial(integratedCoefficients)
}
@ -146,12 +163,11 @@ public fun <C, A> Polynomial<C>.nthAntiderivative(
algebra: A,
order: UInt,
): Polynomial<C> where A : Field<C>, A : NumericAlgebra<C> = algebra {
TODO()
val integratedCoefficients = buildList(coefficients.size + 1) {
add(zero)
coefficients.forEachIndexed{ index, t -> add(t / number(index + 1)) }
val newCoefficients = buildList(coefficients.size + order.toInt()) {
repeat(order.toInt()) { add(zero) }
coefficients.mapIndexedTo(this) { index, c -> (1..order.toInt()).fold(c) { acc, i -> acc / number(index + i) } }
}
Polynomial(integratedCoefficients)
return Polynomial(newCoefficients)
}
/**