Merge branch 'dev' into feature/tensor-algebra

This commit is contained in:
Roland Grinis 2021-03-15 07:29:02 +00:00
commit 7d416f55d4
115 changed files with 2179 additions and 2585 deletions
.github/workflows
CHANGELOG.mdbuild.gradle.kts
examples
gradle.properties
gradle/wrapper
kmath-ast/src
commonMain/kotlin/space/kscience/kmath/ast
jsTest/kotlin/space/kscience/kmath/estree
jvmMain/kotlin/space/kscience/kmath/ast
jvmTest/kotlin/space/kscience/kmath/asm
kmath-commons/src/main/kotlin/space/kscience/kmath/commons
kmath-complex/src
commonMain/kotlin/space/kscience/kmath/complex
commonTest/kotlin/space/kscience/kmath/complex
kmath-core
kmath-coroutines/src/commonMain/kotlin/space/kscience/kmath/chains
kmath-dimensions/src/commonMain/kotlin/space/kscience/kmath/dimensions
kmath-ejml/src/main/kotlin/space/kscience/kmath/ejml
kmath-for-real/src
commonMain/kotlin/space/kscience/kmath/real
commonTest/kotlin/kaceince/kmath/real
kmath-functions/src/commonMain/kotlin/space/kscience/kmath

@ -1,6 +1,5 @@
# Based on https://github.com/touchlab/Stately/blob/main/.github/workflows/build.yml name: Gradle build
name: build
on: [ push ] on: [ push ]
jobs: jobs:

31
.github/workflows/pages.yml vendored Normal file

@ -0,0 +1,31 @@
name: Dokka publication
on:
push:
branches:
- master
jobs:
build:
runs-on: ubuntu-20.04
steps:
- name: Checkout the repo
uses: actions/checkout@v2
- name: Set up JDK 11
uses: actions/setup-java@v1
with:
java-version: 11
- name: Cache gradle
uses: actions/cache@v2
with:
path: ~/.gradle/caches
key: ubuntu-20.04-gradle-${{ hashFiles('*.gradle.kts') }}
restore-keys: |
ubuntu-20.04-gradle-
- name: Build
run: ./gradlew build --no-daemon --stacktrace
- name: Deploy to GitHub Pages
uses: JamesIves/github-pages-deploy-action@4.1.0
with:
branch: gh-pages
folder: build/dokka/htmlMultiModule

@ -1,10 +1,15 @@
# Based on https://github.com/touchlab/Stately/blob/main/.github/workflows/deploy.yml name: Gradle publish
name: deploy on:
on: workflow_dispatch workflow_dispatch:
release:
types:
- created
jobs: jobs:
build: publish:
environment:
name: publish
strategy: strategy:
matrix: matrix:
os: [macOS-latest, windows-latest] os: [macOS-latest, windows-latest]
@ -33,13 +38,22 @@ jobs:
key: ${{ runner.os }}-gradle-${{ hashFiles('*.gradle.kts') }} key: ${{ runner.os }}-gradle-${{ hashFiles('*.gradle.kts') }}
restore-keys: | restore-keys: |
${{ runner.os }}-gradle- ${{ runner.os }}-gradle-
- name: Publish Mac Artifacts
if: matrix.os == 'macOS-latest'
run: ./gradlew publishAllPublicationsToSpaceRepository publishAllPublicationsToSonatypeRepository --no-daemon --stacktrace -PsonatypePublish=true -PsonatypeUser=${{ secrets.SONATYPE_USERNAME }} -PsonatypePassword=${{ secrets.SONATYPE_PASSWORD }} -PspaceUser=${{ secrets.SPACE_USERNAME }} -PspacePassword=${{ secrets.SPACE_PASSWORD }}
env:
signingKey: ${{ secrets.SIGNING_KEY }}
- name: Publish Windows Artifacts - name: Publish Windows Artifacts
if: matrix.os == 'windows-latest' if: matrix.os == 'windows-latest'
run: ./gradlew publishAllPublicationsToSpaceRepository publishAllPublicationsToSonatypeRepository --no-daemon --stacktrace -PsonatypePublish=true -PsonatypeUser=${{ secrets.SONATYPE_USERNAME }} -PsonatypePassword=${{ secrets.SONATYPE_PASSWORD }} -PspaceUser=${{ secrets.SPACE_USERNAME }} -PspacePassword=${{ secrets.SPACE_PASSWORD }} run: >
env: ./gradlew release --no-daemon
signingKey: ${{ secrets.SIGNING_KEY }} -Ppublishing.enabled=true
-Ppublishing.github.user=${{ secrets.PUBLISHING_GITHUB_USER }}
-Ppublishing.github.token=${{ secrets.PUBLISHING_GITHUB_TOKEN }}
-Ppublishing.space.user=${{ secrets.PUBLISHING_SPACE_USER }}
-Ppublishing.space.token=${{ secrets.PUBLISHING_SPACE_TOKEN }}
- name: Publish Mac Artifacts
if: matrix.os == 'macOS-latest'
run: >
./gradlew release --no-daemon
-Ppublishing.enabled=true
-Ppublishing.platform=macosX64
-Ppublishing.github.user=${{ secrets.PUBLISHING_GITHUB_USER }}
-Ppublishing.github.token=${{ secrets.PUBLISHING_GITHUB_TOKEN }}
-Ppublishing.space.user=${{ secrets.PUBLISHING_SPACE_USER }}
-Ppublishing.space.token=${{ secrets.PUBLISHING_SPACE_TOKEN }}

@ -2,12 +2,20 @@
## [Unreleased] ## [Unreleased]
### Added ### Added
- ScaleOperations interface
- Field extends ScaleOperations
### Changed ### Changed
- Exponential operations merged with hyperbolic functions
- Space is replaced by Group. Space is reserved for vector spaces.
- VectorSpace is now a vector space
-
### Deprecated ### Deprecated
### Removed ### Removed
- Nearest in Domain. To be implemented in geometry package.
- Number multiplication and division in main Algebra chain
### Fixed ### Fixed

@ -20,7 +20,7 @@ allprojects {
} }
group = "space.kscience" group = "space.kscience"
version = "0.2.0" version = "0.3.0"
} }
subprojects { subprojects {
@ -32,9 +32,9 @@ readme {
} }
ksciencePublish { ksciencePublish {
spaceRepo = "https://maven.pkg.jetbrains.space/mipt-npm/p/sci/maven" github("kmath")
bintrayRepo = "kscience" space()
githubProject = "kmath" sonatype()
} }
apiValidation { apiValidation {

@ -92,6 +92,14 @@ benchmark {
iterationTimeUnit = "ms" // time unity for iterationTime, default is seconds iterationTimeUnit = "ms" // time unity for iterationTime, default is seconds
include("ExpressionsInterpretersBenchmark") include("ExpressionsInterpretersBenchmark")
} }
configurations.register("matrixInverse") {
warmups = 1 // number of warmup iterations
iterations = 3 // number of iterations
iterationTime = 500 // time in seconds per iteration
iterationTimeUnit = "ms" // time unity for iterationTime, default is seconds
include("MatrixInverseBenchmark")
}
} }
kotlin.sourceSets.all { kotlin.sourceSets.all {

@ -1,34 +1,38 @@
package space.kscience.kmath.benchmarks package space.kscience.kmath.benchmarks
import org.openjdk.jmh.annotations.Benchmark import kotlinx.benchmark.Benchmark
import org.openjdk.jmh.annotations.Scope import kotlinx.benchmark.Blackhole
import org.openjdk.jmh.annotations.State import kotlinx.benchmark.Scope
import kotlinx.benchmark.State
import java.nio.IntBuffer import java.nio.IntBuffer
@State(Scope.Benchmark) @State(Scope.Benchmark)
internal class ArrayBenchmark { internal class ArrayBenchmark {
@Benchmark @Benchmark
fun benchmarkArrayRead() { fun benchmarkArrayRead(blackhole: Blackhole) {
var res = 0 var res = 0
for (i in 1..space.kscience.kmath.benchmarks.ArrayBenchmark.Companion.size) res += space.kscience.kmath.benchmarks.ArrayBenchmark.Companion.array[space.kscience.kmath.benchmarks.ArrayBenchmark.Companion.size - i] for (i in 1..size) res += array[size - i]
blackhole.consume(res)
} }
@Benchmark @Benchmark
fun benchmarkBufferRead() { fun benchmarkBufferRead(blackhole: Blackhole) {
var res = 0 var res = 0
for (i in 1..space.kscience.kmath.benchmarks.ArrayBenchmark.Companion.size) res += space.kscience.kmath.benchmarks.ArrayBenchmark.Companion.arrayBuffer[space.kscience.kmath.benchmarks.ArrayBenchmark.Companion.size - i] for (i in 1..size) res += arrayBuffer[size - i]
blackhole.consume(res)
} }
@Benchmark @Benchmark
fun nativeBufferRead() { fun nativeBufferRead(blackhole: Blackhole) {
var res = 0 var res = 0
for (i in 1..space.kscience.kmath.benchmarks.ArrayBenchmark.Companion.size) res += space.kscience.kmath.benchmarks.ArrayBenchmark.Companion.nativeBuffer[space.kscience.kmath.benchmarks.ArrayBenchmark.Companion.size - i] for (i in 1..size) res += nativeBuffer[size - i]
blackhole.consume(res)
} }
companion object { private companion object {
const val size: Int = 1000 private const val size = 1000
val array: IntArray = IntArray(space.kscience.kmath.benchmarks.ArrayBenchmark.Companion.size) { it } private val array = IntArray(size) { it }
val arrayBuffer: IntBuffer = IntBuffer.wrap(space.kscience.kmath.benchmarks.ArrayBenchmark.Companion.array) private val arrayBuffer = IntBuffer.wrap(array)
val nativeBuffer: IntBuffer = IntBuffer.allocate(space.kscience.kmath.benchmarks.ArrayBenchmark.Companion.size).also { for (i in 0 until space.kscience.kmath.benchmarks.ArrayBenchmark.Companion.size) it.put(i, i) } private val nativeBuffer = IntBuffer.allocate(size).also { for (i in 0 until size) it.put(i, i) }
} }
} }

@ -1,8 +1,8 @@
package space.kscience.kmath.benchmarks package space.kscience.kmath.benchmarks
import org.openjdk.jmh.annotations.Benchmark import kotlinx.benchmark.Benchmark
import org.openjdk.jmh.annotations.Scope import kotlinx.benchmark.Scope
import org.openjdk.jmh.annotations.State import kotlinx.benchmark.State
import space.kscience.kmath.complex.Complex import space.kscience.kmath.complex.Complex
import space.kscience.kmath.complex.complex import space.kscience.kmath.complex.complex
import space.kscience.kmath.structures.MutableBuffer import space.kscience.kmath.structures.MutableBuffer
@ -28,7 +28,7 @@ internal class BufferBenchmark {
} }
} }
companion object { private companion object {
const val size: Int = 100 private const val size = 100
} }
} }

@ -1,67 +1,65 @@
package space.kscience.kmath.benchmarks package space.kscience.kmath.benchmarks
import kotlinx.benchmark.Benchmark import kotlinx.benchmark.Benchmark
import org.openjdk.jmh.annotations.Scope import kotlinx.benchmark.Blackhole
import org.openjdk.jmh.annotations.State import kotlinx.benchmark.Scope
import space.kscience.kmath.commons.linear.CMMatrixContext import kotlinx.benchmark.State
import space.kscience.kmath.ejml.EjmlMatrixContext import space.kscience.kmath.commons.linear.CMLinearSpace
import space.kscience.kmath.linear.BufferMatrixContext import space.kscience.kmath.ejml.EjmlLinearSpace
import space.kscience.kmath.linear.Matrix import space.kscience.kmath.linear.LinearSpace
import space.kscience.kmath.linear.RealMatrixContext import space.kscience.kmath.linear.invoke
import space.kscience.kmath.operations.RealField import space.kscience.kmath.operations.RealField
import space.kscience.kmath.operations.invoke
import space.kscience.kmath.structures.Buffer
import kotlin.random.Random import kotlin.random.Random
@State(Scope.Benchmark) @State(Scope.Benchmark)
internal class DotBenchmark { internal class DotBenchmark {
companion object { companion object {
val random = Random(12224) val random = Random(12224)
val dim = 1000 const val dim = 1000
//creating invertible matrix //creating invertible matrix
val matrix1 = Matrix.real(dim, dim) { i, j -> if (i <= j) random.nextDouble() else 0.0 } val matrix1 = LinearSpace.real.buildMatrix(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 } val matrix2 = LinearSpace.real.buildMatrix(dim, dim) { i, j -> if (i <= j) random.nextDouble() else 0.0 }
val cmMatrix1 = CMMatrixContext { matrix1.toCM() } val cmMatrix1 = CMLinearSpace { matrix1.toCM() }
val cmMatrix2 = CMMatrixContext { matrix2.toCM() } val cmMatrix2 = CMLinearSpace { matrix2.toCM() }
val ejmlMatrix1 = EjmlMatrixContext { matrix1.toEjml() } val ejmlMatrix1 = EjmlLinearSpace { matrix1.toEjml() }
val ejmlMatrix2 = EjmlMatrixContext { matrix2.toEjml() } val ejmlMatrix2 = EjmlLinearSpace { matrix2.toEjml() }
} }
@Benchmark @Benchmark
fun cmDot() { fun cmDot(blackhole: Blackhole) {
CMMatrixContext { CMLinearSpace.run {
cmMatrix1 dot cmMatrix2 blackhole.consume(cmMatrix1 dot cmMatrix2)
} }
} }
@Benchmark @Benchmark
fun ejmlDot() { fun ejmlDot(blackhole: Blackhole) {
EjmlMatrixContext { EjmlLinearSpace {
ejmlMatrix1 dot ejmlMatrix2 blackhole.consume(ejmlMatrix1 dot ejmlMatrix2)
} }
} }
@Benchmark @Benchmark
fun ejmlDotWithConversion() { fun ejmlDotWithConversion(blackhole: Blackhole) {
EjmlMatrixContext { EjmlLinearSpace {
matrix1 dot matrix2 blackhole.consume(matrix1 dot matrix2)
} }
} }
@Benchmark @Benchmark
fun bufferedDot() { fun bufferedDot(blackhole: Blackhole) {
BufferMatrixContext(RealField, Buffer.Companion::real).invoke { LinearSpace.auto(RealField).invoke {
matrix1 dot matrix2 blackhole.consume(matrix1 dot matrix2)
} }
} }
@Benchmark @Benchmark
fun realDot() { fun realDot(blackhole: Blackhole) {
RealMatrixContext { LinearSpace.real {
matrix1 dot matrix2 blackhole.consume(matrix1 dot matrix2)
} }
} }
} }

@ -1,64 +1,62 @@
package space.kscience.kmath.benchmarks package space.kscience.kmath.benchmarks
import org.openjdk.jmh.annotations.Benchmark import kotlinx.benchmark.Benchmark
import org.openjdk.jmh.annotations.Scope import kotlinx.benchmark.Blackhole
import org.openjdk.jmh.annotations.State import kotlinx.benchmark.Scope
import kotlinx.benchmark.State
import space.kscience.kmath.asm.compile import space.kscience.kmath.asm.compile
import space.kscience.kmath.ast.mstInField import space.kscience.kmath.ast.mstInField
import space.kscience.kmath.expressions.Expression import space.kscience.kmath.expressions.Expression
import space.kscience.kmath.expressions.expressionInField import space.kscience.kmath.expressions.expressionInField
import space.kscience.kmath.expressions.invoke import space.kscience.kmath.expressions.invoke
import space.kscience.kmath.expressions.symbol import space.kscience.kmath.expressions.symbol
import space.kscience.kmath.operations.Field
import space.kscience.kmath.operations.RealField import space.kscience.kmath.operations.RealField
import space.kscience.kmath.operations.bindSymbol import space.kscience.kmath.operations.bindSymbol
import kotlin.random.Random import kotlin.random.Random
@State(Scope.Benchmark) @State(Scope.Benchmark)
internal class ExpressionsInterpretersBenchmark { internal class ExpressionsInterpretersBenchmark {
private val algebra: Field<Double> = RealField
val x by symbol
@Benchmark @Benchmark
fun functionalExpression() { fun functionalExpression(blackhole: Blackhole) {
val expr = algebra.expressionInField { val expr = algebra.expressionInField {
val x = bindSymbol(x) val x = bindSymbol(x)
x * const(2.0) + const(2.0) / x - const(16.0) x * const(2.0) + const(2.0) / x - const(16.0)
} }
invokeAndSum(expr) invokeAndSum(expr, blackhole)
} }
@Benchmark @Benchmark
fun mstExpression() { fun mstExpression(blackhole: Blackhole) {
val expr = algebra.mstInField { val expr = algebra.mstInField {
val x = bindSymbol(x) val x = bindSymbol(x)
x * 2.0 + 2.0 / x - 16.0 x * 2.0 + number(2.0) / x - 16.0
} }
invokeAndSum(expr) invokeAndSum(expr, blackhole)
} }
@Benchmark @Benchmark
fun asmExpression() { fun asmExpression(blackhole: Blackhole) {
val expr = algebra.mstInField { val expr = algebra.mstInField {
val x = bindSymbol(x) val x = bindSymbol(x)
x * 2.0 + 2.0 / x - 16.0 x * 2.0 + number(2.0) / x - 16.0
}.compile() }.compile()
invokeAndSum(expr) invokeAndSum(expr, blackhole)
} }
@Benchmark @Benchmark
fun rawExpression() { fun rawExpression(blackhole: Blackhole) {
val expr = Expression<Double> { args -> val expr = Expression<Double> { args ->
val x = args.getValue(x) val x = args.getValue(x)
x * 2.0 + 2.0 / x - 16.0 x * 2.0 + 2.0 / x - 16.0
} }
invokeAndSum(expr)
invokeAndSum(expr, blackhole)
} }
private fun invokeAndSum(expr: Expression<Double>) { private fun invokeAndSum(expr: Expression<Double>, blackhole: Blackhole) {
val random = Random(0) val random = Random(0)
var sum = 0.0 var sum = 0.0
@ -66,6 +64,11 @@ internal class ExpressionsInterpretersBenchmark {
sum += expr(x to random.nextDouble()) sum += expr(x to random.nextDouble())
} }
println(sum) blackhole.consume(sum)
}
private companion object {
private val algebra = RealField
private val x by symbol
} }
} }

@ -1,47 +0,0 @@
package space.kscience.kmath.benchmarks
import kotlinx.benchmark.Benchmark
import org.openjdk.jmh.annotations.Scope
import org.openjdk.jmh.annotations.State
import space.kscience.kmath.commons.linear.CMMatrixContext
import space.kscience.kmath.commons.linear.CMMatrixContext.dot
import space.kscience.kmath.commons.linear.inverse
import space.kscience.kmath.ejml.EjmlMatrixContext
import space.kscience.kmath.ejml.inverse
import space.kscience.kmath.linear.Matrix
import space.kscience.kmath.linear.MatrixContext
import space.kscience.kmath.linear.inverseWithLup
import space.kscience.kmath.linear.real
import kotlin.random.Random
@State(Scope.Benchmark)
internal class LinearAlgebraBenchmark {
companion object {
val random = Random(1224)
val dim = 100
//creating invertible matrix
val u = Matrix.real(dim, dim) { i, j -> if (i <= j) random.nextDouble() else 0.0 }
val l = Matrix.real(dim, dim) { i, j -> if (i >= j) random.nextDouble() else 0.0 }
val matrix = l dot u
}
@Benchmark
fun kmathLupInversion() {
MatrixContext.real.inverseWithLup(matrix)
}
@Benchmark
fun cmLUPInversion() {
with(CMMatrixContext) {
inverse(matrix)
}
}
@Benchmark
fun ejmlInverse() {
with(EjmlMatrixContext) {
inverse(matrix)
}
}
}

@ -0,0 +1,48 @@
package space.kscience.kmath.benchmarks
import kotlinx.benchmark.Benchmark
import kotlinx.benchmark.Blackhole
import kotlinx.benchmark.Scope
import kotlinx.benchmark.State
import space.kscience.kmath.commons.linear.CMLinearSpace
import space.kscience.kmath.commons.linear.inverse
import space.kscience.kmath.ejml.EjmlLinearSpace
import space.kscience.kmath.ejml.inverse
import space.kscience.kmath.linear.LinearSpace
import space.kscience.kmath.linear.inverseWithLup
import space.kscience.kmath.linear.invoke
import kotlin.random.Random
@State(Scope.Benchmark)
internal class MatrixInverseBenchmark {
companion object {
val random = Random(1224)
const val dim = 100
private val space = LinearSpace.real
//creating invertible matrix
val u = space.buildMatrix(dim, dim) { i, j -> if (i <= j) random.nextDouble() else 0.0 }
val l = space.buildMatrix(dim, dim) { i, j -> if (i >= j) random.nextDouble() else 0.0 }
val matrix = space { l dot u }
}
@Benchmark
fun kmathLupInversion(blackhole: Blackhole) {
blackhole.consume(LinearSpace.real.inverseWithLup(matrix))
}
@Benchmark
fun cmLUPInversion(blackhole: Blackhole) {
with(CMLinearSpace) {
blackhole.consume(inverse(matrix))
}
}
@Benchmark
fun ejmlInverse(blackhole: Blackhole) {
with(EjmlLinearSpace) {
blackhole.consume(inverse(matrix))
}
}
}

@ -1,8 +1,9 @@
package space.kscience.kmath.benchmarks package space.kscience.kmath.benchmarks
import org.openjdk.jmh.annotations.Benchmark import kotlinx.benchmark.Benchmark
import org.openjdk.jmh.annotations.Scope import kotlinx.benchmark.Blackhole
import org.openjdk.jmh.annotations.State import kotlinx.benchmark.Scope
import kotlinx.benchmark.State
import space.kscience.kmath.nd.* import space.kscience.kmath.nd.*
import space.kscience.kmath.operations.RealField import space.kscience.kmath.operations.RealField
import space.kscience.kmath.structures.Buffer import space.kscience.kmath.structures.Buffer
@ -10,35 +11,38 @@ import space.kscience.kmath.structures.Buffer
@State(Scope.Benchmark) @State(Scope.Benchmark)
internal class NDFieldBenchmark { internal class NDFieldBenchmark {
@Benchmark @Benchmark
fun autoFieldAdd() { fun autoFieldAdd(blackhole: Blackhole) {
with(autoField) { with(autoField) {
var res: NDStructure<Double> = one var res: NDStructure<Double> = one
repeat(n) { res += one } repeat(n) { res += one }
blackhole.consume(res)
} }
} }
@Benchmark @Benchmark
fun specializedFieldAdd() { fun specializedFieldAdd(blackhole: Blackhole) {
with(specializedField) { with(specializedField) {
var res: NDStructure<Double> = one var res: NDStructure<Double> = one
repeat(n) { res += 1.0 } repeat(n) { res += 1.0 }
blackhole.consume(res)
} }
} }
@Benchmark @Benchmark
fun boxingFieldAdd() { fun boxingFieldAdd(blackhole: Blackhole) {
with(genericField) { with(genericField) {
var res: NDStructure<Double> = one var res: NDStructure<Double> = one
repeat(n) { res += 1.0 } repeat(n) { res += 1.0 }
blackhole.consume(res)
} }
} }
companion object { private companion object {
const val dim: Int = 1000 private const val dim = 1000
const val n: Int = 100 private const val n = 100
val autoField = NDAlgebra.auto(RealField, dim, dim) private val autoField = NDAlgebra.auto(RealField, dim, dim)
val specializedField: RealNDField = NDAlgebra.real(dim, dim) private val specializedField = NDAlgebra.real(dim, dim)
val genericField = NDAlgebra.field(RealField, Buffer.Companion::boxing, dim, dim) private val genericField = NDAlgebra.field(RealField, Buffer.Companion::boxing, dim, dim)
} }
} }

@ -1,51 +1,61 @@
package space.kscience.kmath.benchmarks package space.kscience.kmath.benchmarks
import kotlinx.benchmark.Benchmark
import kotlinx.benchmark.Blackhole
import kotlinx.benchmark.Scope
import kotlinx.benchmark.State
import org.jetbrains.bio.viktor.F64Array import org.jetbrains.bio.viktor.F64Array
import org.openjdk.jmh.annotations.Benchmark import space.kscience.kmath.nd.NDAlgebra
import org.openjdk.jmh.annotations.Scope import space.kscience.kmath.nd.NDStructure
import org.openjdk.jmh.annotations.State import space.kscience.kmath.nd.auto
import space.kscience.kmath.nd.* import space.kscience.kmath.nd.real
import space.kscience.kmath.operations.RealField import space.kscience.kmath.operations.RealField
import space.kscience.kmath.viktor.ViktorNDField import space.kscience.kmath.viktor.ViktorNDField
@State(Scope.Benchmark) @State(Scope.Benchmark)
internal class ViktorBenchmark { internal class ViktorBenchmark {
final val dim: Int = 1000
final val n: Int = 100
// automatically build context most suited for given type.
final val autoField: NDField<Double, RealField> = NDAlgebra.auto(RealField, dim, dim)
final val realField: RealNDField = NDAlgebra.real(dim, dim)
final val viktorField: ViktorNDField = ViktorNDField(dim, dim)
@Benchmark @Benchmark
fun automaticFieldAddition() { fun automaticFieldAddition(blackhole: Blackhole) {
with(autoField) { with(autoField) {
var res: NDStructure<Double> = one var res: NDStructure<Double> = one
repeat(n) { res += 1.0 } repeat(n) { res += 1.0 }
blackhole.consume(res)
} }
} }
@Benchmark @Benchmark
fun realFieldAddition() { fun realFieldAddition(blackhole: Blackhole) {
with(realField) { with(realField) {
var res: NDStructure<Double> = one var res: NDStructure<Double> = one
repeat(n) { res += 1.0 } repeat(n) { res += 1.0 }
blackhole.consume(res)
} }
} }
@Benchmark @Benchmark
fun viktorFieldAddition() { fun viktorFieldAddition(blackhole: Blackhole) {
with(viktorField) { with(viktorField) {
var res = one var res = one
repeat(n) { res += 1.0 } repeat(n) { res += 1.0 }
blackhole.consume(res)
} }
} }
@Benchmark @Benchmark
fun rawViktor() { fun rawViktor(blackhole: Blackhole) {
val one = F64Array.full(init = 1.0, shape = intArrayOf(dim, dim)) val one = F64Array.full(init = 1.0, shape = intArrayOf(dim, dim))
var res = one var res = one
repeat(n) { res = res + one } repeat(n) { res = res + one }
blackhole.consume(res)
}
private companion object {
private const val dim = 1000
private const val n = 100
// automatically build context most suited for given type.
private val autoField = NDAlgebra.auto(RealField, dim, dim)
private val realField = NDAlgebra.real(dim, dim)
private val viktorField = ViktorNDField(dim, dim)
} }
} }

@ -1,48 +1,53 @@
package space.kscience.kmath.benchmarks package space.kscience.kmath.benchmarks
import kotlinx.benchmark.Benchmark
import kotlinx.benchmark.Blackhole
import kotlinx.benchmark.Scope
import kotlinx.benchmark.State
import org.jetbrains.bio.viktor.F64Array import org.jetbrains.bio.viktor.F64Array
import org.openjdk.jmh.annotations.Benchmark import space.kscience.kmath.nd.NDAlgebra
import org.openjdk.jmh.annotations.Scope import space.kscience.kmath.nd.auto
import org.openjdk.jmh.annotations.State import space.kscience.kmath.nd.real
import space.kscience.kmath.nd.*
import space.kscience.kmath.operations.RealField import space.kscience.kmath.operations.RealField
import space.kscience.kmath.viktor.ViktorNDField import space.kscience.kmath.viktor.ViktorNDField
@State(Scope.Benchmark) @State(Scope.Benchmark)
internal class ViktorLogBenchmark { internal class ViktorLogBenchmark {
final val dim: Int = 1000
final val n: Int = 100
// automatically build context most suited for given type.
final val autoField: NDField<Double, RealField> = NDAlgebra.auto(RealField, dim, dim)
final val realField: RealNDField = NDAlgebra.real(dim, dim)
final val viktorField: ViktorNDField = ViktorNDField(intArrayOf(dim, dim))
@Benchmark @Benchmark
fun realFieldLog() { fun realFieldLog(blackhole: Blackhole) {
with(realField) { with(realNdField) {
val fortyTwo = produce { 42.0 } val fortyTwo = produce { 42.0 }
var res = one var res = one
repeat(n) { res = ln(fortyTwo) } repeat(n) { res = ln(fortyTwo) }
blackhole.consume(res)
} }
} }
@Benchmark @Benchmark
fun viktorFieldLog() { fun viktorFieldLog(blackhole: Blackhole) {
with(viktorField) { with(viktorField) {
val fortyTwo = produce { 42.0 } val fortyTwo = produce { 42.0 }
var res = one var res = one
repeat(n) { res = ln(fortyTwo) } repeat(n) { res = ln(fortyTwo) }
blackhole.consume(res)
} }
} }
@Benchmark @Benchmark
fun rawViktorLog() { fun rawViktorLog(blackhole: Blackhole) {
val fortyTwo = F64Array.full(dim, dim, init = 42.0) val fortyTwo = F64Array.full(dim, dim, init = 42.0)
var res: F64Array lateinit var res: F64Array
repeat(n) { repeat(n) { res = fortyTwo.log() }
res = fortyTwo.log() blackhole.consume(res)
} }
private companion object {
private const val dim = 1000
private const val n = 100
// automatically build context most suited for given type.
private val autoField = NDAlgebra.auto(RealField, dim, dim)
private val realNdField = NDAlgebra.real(dim, dim)
private val viktorField = ViktorNDField(intArrayOf(dim, dim))
} }
} }

@ -6,10 +6,10 @@ import space.kscience.kmath.operations.RealField
fun main() { fun main() {
val expr = RealField.mstInField { val expr = RealField.mstInField {
val x = bindSymbol("x") val x = bindSymbol("x")
x * 2.0 + 2.0 / x - 16.0 x * 2.0 + number(2.0) / x - 16.0
} }
repeat(10000000){ repeat(10000000) {
expr.invoke("x" to 1.0) expr.invoke("x" to 1.0)
} }
} }

@ -0,0 +1,28 @@
package space.kscience.kmath.linear
import space.kscience.kmath.real.*
import space.kscience.kmath.structures.RealBuffer
fun main() {
val x0 = Point(0.0, 0.0, 0.0)
val sigma = Point(1.0, 1.0, 1.0)
val gaussian: (Point<Double>) -> Double = { x ->
require(x.size == x0.size)
kotlin.math.exp(-((x - x0) / sigma).square().sum())
}
fun ((Point<Double>) -> Double).grad(x: Point<Double>): Point<Double> {
require(x.size == x0.size)
return RealBuffer(x.size) { i ->
val h = sigma[i] / 5
val dVector = RealBuffer(x.size) { if (it == i) h else 0.0 }
val f1 = invoke(x + dVector / 2)
val f0 = invoke(x - dVector / 2)
(f1 - f0) / h
}
}
println(gaussian.grad(x0))
}

@ -3,8 +3,8 @@ package space.kscience.kmath.structures
import space.kscience.kmath.misc.UnstableKMathAPI import space.kscience.kmath.misc.UnstableKMathAPI
import space.kscience.kmath.nd.* import space.kscience.kmath.nd.*
import space.kscience.kmath.operations.ExtendedField import space.kscience.kmath.operations.ExtendedField
import space.kscience.kmath.operations.NumbersAddOperations
import space.kscience.kmath.operations.RealField import space.kscience.kmath.operations.RealField
import space.kscience.kmath.operations.RingWithNumbers
import java.util.* import java.util.*
import java.util.stream.IntStream import java.util.stream.IntStream
@ -15,7 +15,7 @@ import java.util.stream.IntStream
class StreamRealNDField( class StreamRealNDField(
override val shape: IntArray, override val shape: IntArray,
) : NDField<Double, RealField>, ) : NDField<Double, RealField>,
RingWithNumbers<NDStructure<Double>>, NumbersAddOperations<NDStructure<Double>>,
ExtendedField<NDStructure<Double>> { ExtendedField<NDStructure<Double>> {
private val strides = DefaultStrides(shape) private val strides = DefaultStrides(shape)
@ -79,25 +79,29 @@ class StreamRealNDField(
return NDBuffer(strides, array.asBuffer()) return NDBuffer(strides, array.asBuffer())
} }
override fun power(arg: NDStructure<Double>, pow: Number): NDBuffer<Double> = arg.map() { power(it, pow) } override fun NDStructure<Double>.unaryMinus(): NDStructure<Double> = map { -it }
override fun exp(arg: NDStructure<Double>): NDBuffer<Double> = arg.map() { exp(it) } override fun scale(a: NDStructure<Double>, value: Double): NDStructure<Double> = a.map { it * value }
override fun ln(arg: NDStructure<Double>): NDBuffer<Double> = arg.map() { ln(it) } override fun power(arg: NDStructure<Double>, pow: Number): NDBuffer<Double> = arg.map { power(it, pow) }
override fun sin(arg: NDStructure<Double>): NDBuffer<Double> = arg.map() { sin(it) } override fun exp(arg: NDStructure<Double>): NDBuffer<Double> = arg.map { exp(it) }
override fun cos(arg: NDStructure<Double>): NDBuffer<Double> = arg.map() { cos(it) }
override fun tan(arg: NDStructure<Double>): NDBuffer<Double> = arg.map() { tan(it) }
override fun asin(arg: NDStructure<Double>): NDBuffer<Double> = arg.map() { asin(it) }
override fun acos(arg: NDStructure<Double>): NDBuffer<Double> = arg.map() { acos(it) }
override fun atan(arg: NDStructure<Double>): NDBuffer<Double> = arg.map() { atan(it) }
override fun sinh(arg: NDStructure<Double>): NDBuffer<Double> = arg.map() { sinh(it) } override fun ln(arg: NDStructure<Double>): NDBuffer<Double> = arg.map { ln(it) }
override fun cosh(arg: NDStructure<Double>): NDBuffer<Double> = arg.map() { cosh(it) }
override fun tanh(arg: NDStructure<Double>): NDBuffer<Double> = arg.map() { tanh(it) } override fun sin(arg: NDStructure<Double>): NDBuffer<Double> = arg.map { sin(it) }
override fun asinh(arg: NDStructure<Double>): NDBuffer<Double> = arg.map() { asinh(it) } override fun cos(arg: NDStructure<Double>): NDBuffer<Double> = arg.map { cos(it) }
override fun acosh(arg: NDStructure<Double>): NDBuffer<Double> = arg.map() { acosh(it) } override fun tan(arg: NDStructure<Double>): NDBuffer<Double> = arg.map { tan(it) }
override fun atanh(arg: NDStructure<Double>): NDBuffer<Double> = arg.map() { atanh(it) } override fun asin(arg: NDStructure<Double>): NDBuffer<Double> = arg.map { asin(it) }
override fun acos(arg: NDStructure<Double>): NDBuffer<Double> = arg.map { acos(it) }
override fun atan(arg: NDStructure<Double>): NDBuffer<Double> = arg.map { atan(it) }
override fun sinh(arg: NDStructure<Double>): NDBuffer<Double> = arg.map { sinh(it) }
override fun cosh(arg: NDStructure<Double>): NDBuffer<Double> = arg.map { cosh(it) }
override fun tanh(arg: NDStructure<Double>): NDBuffer<Double> = arg.map { tanh(it) }
override fun asinh(arg: NDStructure<Double>): NDBuffer<Double> = arg.map { asinh(it) }
override fun acosh(arg: NDStructure<Double>): NDBuffer<Double> = arg.map { acosh(it) }
override fun atanh(arg: NDStructure<Double>): NDBuffer<Double> = arg.map { atanh(it) }
} }
fun NDAlgebra.Companion.realWithStream(vararg shape: Int): StreamRealNDField = StreamRealNDField(shape) fun NDAlgebra.Companion.realWithStream(vararg shape: Int): StreamRealNDField = StreamRealNDField(shape)

@ -7,7 +7,7 @@ import kotlin.system.measureTimeMillis
@Suppress("UNUSED_VARIABLE") @Suppress("UNUSED_VARIABLE")
fun main() { fun main() {
val n = 6000 val n = 6000
val structure = NDStructure.build(intArrayOf(n, n), Buffer.Companion::auto) { 1.0 } val structure = NDStructure.buffered(intArrayOf(n, n), Buffer.Companion::auto) { 1.0 }
structure.mapToBuffer { it + 1 } // warm-up structure.mapToBuffer { it + 1 } // warm-up
val time1 = measureTimeMillis { val res = structure.mapToBuffer { it + 1 } } val time1 = measureTimeMillis { val res = structure.mapToBuffer { it + 1 } }
println("Structure mapping finished in $time1 millis") println("Structure mapping finished in $time1 millis")

@ -5,7 +5,7 @@ import space.kscience.kmath.dimensions.D3
import space.kscience.kmath.dimensions.DMatrixContext import space.kscience.kmath.dimensions.DMatrixContext
import space.kscience.kmath.dimensions.Dimension import space.kscience.kmath.dimensions.Dimension
private fun DMatrixContext<Double>.simple() { private fun DMatrixContext<Double, *>.simple() {
val m1 = produce<D2, D3> { i, j -> (i + j).toDouble() } val m1 = produce<D2, D3> { i, j -> (i + j).toDouble() }
val m2 = produce<D3, D2> { i, j -> (i + j).toDouble() } val m2 = produce<D3, D2> { i, j -> (i + j).toDouble() }
@ -17,7 +17,7 @@ private object D5 : Dimension {
override val dim: UInt = 5u override val dim: UInt = 5u
} }
private fun DMatrixContext<Double>.custom() { private fun DMatrixContext<Double, *>.custom() {
val m1 = produce<D2, D5> { i, j -> (i + j).toDouble() } val m1 = produce<D2, D5> { i, j -> (i + j).toDouble() }
val m2 = produce<D5, D2> { i, j -> (i - j).toDouble() } val m2 = produce<D5, D2> { i, j -> (i - j).toDouble() }
val m3 = produce<D2, D2> { i, j -> (i - j).toDouble() } val m3 = produce<D2, D2> { i, j -> (i - j).toDouble() }

@ -1,9 +1,10 @@
kotlin.code.style=official kotlin.code.style=official
kotlin.mpp.enableGranularSourceSetsMetadata=true
kotlin.mpp.stability.nowarn=true kotlin.mpp.stability.nowarn=true
kotlin.native.enableDependencyPropagation=false
kotlin.parallel.tasks.in.project=true kotlin.parallel.tasks.in.project=true
org.gradle.configureondemand=true org.gradle.configureondemand=true
org.gradle.jvmargs=-XX:MaxMetaspaceSize=512m org.gradle.jvmargs=-XX:MaxMetaspaceSize=1G
org.gradle.parallel=true org.gradle.parallel=true
systemProp.org.gradle.internal.publish.checksums.insecure=true
kotlin.mpp.enableGranularSourceSetsMetadata=true
kotlin.native.enableDependencyPropagation=false

@ -1,5 +1,5 @@
distributionBase=GRADLE_USER_HOME distributionBase=GRADLE_USER_HOME
distributionPath=wrapper/dists distributionPath=wrapper/dists
distributionUrl=https\://services.gradle.org/distributions/gradle-6.8.2-bin.zip distributionUrl=https\://services.gradle.org/distributions/gradle-6.8.3-bin.zip
zipStoreBase=GRADLE_USER_HOME zipStoreBase=GRADLE_USER_HOME
zipStorePath=wrapper/dists zipStorePath=wrapper/dists

@ -18,25 +18,25 @@ public object MstAlgebra : NumericAlgebra<MST> {
} }
/** /**
* [Space] over [MST] nodes. * [Group] over [MST] nodes.
*/ */
public object MstSpace : Space<MST>, NumericAlgebra<MST> { public object MstGroup : Group<MST>, NumericAlgebra<MST>, ScaleOperations<MST> {
public override val zero: MST.Numeric by lazy { number(0.0) } public override val zero: MST.Numeric = number(0.0)
public override fun number(value: Number): MST.Numeric = MstAlgebra.number(value) public override fun number(value: Number): MST.Numeric = MstAlgebra.number(value)
public override fun bindSymbol(value: String): MST.Symbolic = MstAlgebra.bindSymbol(value) public override fun bindSymbol(value: String): MST.Symbolic = MstAlgebra.bindSymbol(value)
public override fun add(a: MST, b: MST): MST.Binary = binaryOperationFunction(SpaceOperations.PLUS_OPERATION)(a, b) public override fun add(a: MST, b: MST): MST.Binary = binaryOperationFunction(GroupOperations.PLUS_OPERATION)(a, b)
public override operator fun MST.unaryPlus(): MST.Unary = public override operator fun MST.unaryPlus(): MST.Unary =
unaryOperationFunction(SpaceOperations.PLUS_OPERATION)(this) unaryOperationFunction(GroupOperations.PLUS_OPERATION)(this)
public override operator fun MST.unaryMinus(): MST.Unary = public override operator fun MST.unaryMinus(): MST.Unary =
unaryOperationFunction(SpaceOperations.MINUS_OPERATION)(this) unaryOperationFunction(GroupOperations.MINUS_OPERATION)(this)
public override operator fun MST.minus(b: MST): MST.Binary = public override operator fun MST.minus(b: MST): MST.Binary =
binaryOperationFunction(SpaceOperations.MINUS_OPERATION)(this, b) binaryOperationFunction(GroupOperations.MINUS_OPERATION)(this, b)
public override fun multiply(a: MST, k: Number): MST.Binary = public override fun scale(a: MST, value: Double): MST.Binary =
binaryOperationFunction(RingOperations.TIMES_OPERATION)(a, number(k)) binaryOperationFunction(RingOperations.TIMES_OPERATION)(a, number(value))
public override fun binaryOperationFunction(operation: String): (left: MST, right: MST) -> MST.Binary = public override fun binaryOperationFunction(operation: String): (left: MST, right: MST) -> MST.Binary =
MstAlgebra.binaryOperationFunction(operation) MstAlgebra.binaryOperationFunction(operation)
@ -49,25 +49,26 @@ public object MstSpace : Space<MST>, NumericAlgebra<MST> {
* [Ring] over [MST] nodes. * [Ring] over [MST] nodes.
*/ */
@OptIn(UnstableKMathAPI::class) @OptIn(UnstableKMathAPI::class)
public object MstRing : Ring<MST>, RingWithNumbers<MST> { public object MstRing : Ring<MST>, NumbersAddOperations<MST>, ScaleOperations<MST> {
public override val zero: MST.Numeric public override val zero: MST.Numeric get() = MstGroup.zero
get() = MstSpace.zero public override val one: MST.Numeric = number(1.0)
public override val one: MST.Numeric by lazy { number(1.0) } public override fun number(value: Number): MST.Numeric = MstGroup.number(value)
public override fun bindSymbol(value: String): MST.Symbolic = MstAlgebra.bindSymbol(value)
public override fun add(a: MST, b: MST): MST.Binary = MstGroup.add(a, b)
public override fun scale(a: MST, value: Double): MST.Binary =
MstGroup.binaryOperationFunction(RingOperations.TIMES_OPERATION)(a, MstGroup.number(value))
public override fun number(value: Number): MST.Numeric = MstSpace.number(value)
public override fun bindSymbol(value: String): MST.Symbolic = MstSpace.bindSymbol(value)
public override fun add(a: MST, b: MST): MST.Binary = MstSpace.add(a, b)
public override fun multiply(a: MST, k: Number): MST.Binary = MstSpace.multiply(a, k)
public override fun multiply(a: MST, b: MST): MST.Binary = public override fun multiply(a: MST, b: MST): MST.Binary =
binaryOperationFunction(RingOperations.TIMES_OPERATION)(a, b) binaryOperationFunction(RingOperations.TIMES_OPERATION)(a, b)
public override operator fun MST.unaryPlus(): MST.Unary = MstSpace { +this@unaryPlus } public override operator fun MST.unaryPlus(): MST.Unary = MstGroup { +this@unaryPlus }
public override operator fun MST.unaryMinus(): MST.Unary = MstSpace { -this@unaryMinus } public override operator fun MST.unaryMinus(): MST.Unary = MstGroup { -this@unaryMinus }
public override operator fun MST.minus(b: MST): MST.Binary = MstSpace { this@minus - b } public override operator fun MST.minus(b: MST): MST.Binary = MstGroup { this@minus - b }
public override fun binaryOperationFunction(operation: String): (left: MST, right: MST) -> MST.Binary = public override fun binaryOperationFunction(operation: String): (left: MST, right: MST) -> MST.Binary =
MstSpace.binaryOperationFunction(operation) MstGroup.binaryOperationFunction(operation)
public override fun unaryOperationFunction(operation: String): (arg: MST) -> MST.Unary = public override fun unaryOperationFunction(operation: String): (arg: MST) -> MST.Unary =
MstAlgebra.unaryOperationFunction(operation) MstAlgebra.unaryOperationFunction(operation)
@ -77,17 +78,18 @@ public object MstRing : Ring<MST>, RingWithNumbers<MST> {
* [Field] over [MST] nodes. * [Field] over [MST] nodes.
*/ */
@OptIn(UnstableKMathAPI::class) @OptIn(UnstableKMathAPI::class)
public object MstField : Field<MST>, RingWithNumbers<MST> { public object MstField : Field<MST>, NumbersAddOperations<MST>, ScaleOperations<MST> {
public override val zero: MST.Numeric public override val zero: MST.Numeric get() = MstRing.zero
get() = MstRing.zero
public override val one: MST.Numeric public override val one: MST.Numeric get() = MstRing.one
get() = MstRing.one
public override fun bindSymbol(value: String): MST.Symbolic = MstRing.bindSymbol(value) public override fun bindSymbol(value: String): MST.Symbolic = MstAlgebra.bindSymbol(value)
public override fun number(value: Number): MST.Numeric = MstRing.number(value) public override fun number(value: Number): MST.Numeric = MstRing.number(value)
public override fun add(a: MST, b: MST): MST.Binary = MstRing.add(a, b) public override fun add(a: MST, b: MST): MST.Binary = MstRing.add(a, b)
public override fun multiply(a: MST, k: Number): MST.Binary = MstRing.multiply(a, k)
public override fun scale(a: MST, value: Double): MST.Binary =
MstGroup.binaryOperationFunction(RingOperations.TIMES_OPERATION)(a, MstGroup.number(value))
public override fun multiply(a: MST, b: MST): MST.Binary = MstRing.multiply(a, b) public override fun multiply(a: MST, b: MST): MST.Binary = MstRing.multiply(a, b)
public override fun divide(a: MST, b: MST): MST.Binary = public override fun divide(a: MST, b: MST): MST.Binary =
binaryOperationFunction(FieldOperations.DIV_OPERATION)(a, b) binaryOperationFunction(FieldOperations.DIV_OPERATION)(a, b)
@ -107,13 +109,10 @@ public object MstField : Field<MST>, RingWithNumbers<MST> {
* [ExtendedField] over [MST] nodes. * [ExtendedField] over [MST] nodes.
*/ */
public object MstExtendedField : ExtendedField<MST>, NumericAlgebra<MST> { public object MstExtendedField : ExtendedField<MST>, NumericAlgebra<MST> {
public override val zero: MST.Numeric public override val zero: MST.Numeric get() = MstField.zero
get() = MstField.zero public override val one: MST.Numeric get() = MstField.one
public override val one: MST.Numeric public override fun bindSymbol(value: String): MST.Symbolic = MstAlgebra.bindSymbol(value)
get() = MstField.one
public override fun bindSymbol(value: String): MST.Symbolic = MstField.bindSymbol(value)
public override fun number(value: Number): MST.Numeric = MstRing.number(value) public override fun number(value: Number): MST.Numeric = MstRing.number(value)
public override fun sin(arg: MST): MST.Unary = unaryOperationFunction(TrigonometricOperations.SIN_OPERATION)(arg) public override fun sin(arg: MST): MST.Unary = unaryOperationFunction(TrigonometricOperations.SIN_OPERATION)(arg)
public override fun cos(arg: MST): MST.Unary = unaryOperationFunction(TrigonometricOperations.COS_OPERATION)(arg) public override fun cos(arg: MST): MST.Unary = unaryOperationFunction(TrigonometricOperations.COS_OPERATION)(arg)
@ -121,14 +120,17 @@ public object MstExtendedField : ExtendedField<MST>, NumericAlgebra<MST> {
public override fun asin(arg: MST): MST.Unary = unaryOperationFunction(TrigonometricOperations.ASIN_OPERATION)(arg) public override fun asin(arg: MST): MST.Unary = unaryOperationFunction(TrigonometricOperations.ASIN_OPERATION)(arg)
public override fun acos(arg: MST): MST.Unary = unaryOperationFunction(TrigonometricOperations.ACOS_OPERATION)(arg) public override fun acos(arg: MST): MST.Unary = unaryOperationFunction(TrigonometricOperations.ACOS_OPERATION)(arg)
public override fun atan(arg: MST): MST.Unary = unaryOperationFunction(TrigonometricOperations.ATAN_OPERATION)(arg) public override fun atan(arg: MST): MST.Unary = unaryOperationFunction(TrigonometricOperations.ATAN_OPERATION)(arg)
public override fun sinh(arg: MST): MST.Unary = unaryOperationFunction(HyperbolicOperations.SINH_OPERATION)(arg) public override fun sinh(arg: MST): MST.Unary = unaryOperationFunction(ExponentialOperations.SINH_OPERATION)(arg)
public override fun cosh(arg: MST): MST.Unary = unaryOperationFunction(HyperbolicOperations.COSH_OPERATION)(arg) public override fun cosh(arg: MST): MST.Unary = unaryOperationFunction(ExponentialOperations.COSH_OPERATION)(arg)
public override fun tanh(arg: MST): MST.Unary = unaryOperationFunction(HyperbolicOperations.TANH_OPERATION)(arg) public override fun tanh(arg: MST): MST.Unary = unaryOperationFunction(ExponentialOperations.TANH_OPERATION)(arg)
public override fun asinh(arg: MST): MST.Unary = unaryOperationFunction(HyperbolicOperations.ASINH_OPERATION)(arg) public override fun asinh(arg: MST): MST.Unary = unaryOperationFunction(ExponentialOperations.ASINH_OPERATION)(arg)
public override fun acosh(arg: MST): MST.Unary = unaryOperationFunction(HyperbolicOperations.ACOSH_OPERATION)(arg) public override fun acosh(arg: MST): MST.Unary = unaryOperationFunction(ExponentialOperations.ACOSH_OPERATION)(arg)
public override fun atanh(arg: MST): MST.Unary = unaryOperationFunction(HyperbolicOperations.ATANH_OPERATION)(arg) public override fun atanh(arg: MST): MST.Unary = unaryOperationFunction(ExponentialOperations.ATANH_OPERATION)(arg)
public override fun add(a: MST, b: MST): MST.Binary = MstField.add(a, b) public override fun add(a: MST, b: MST): MST.Binary = MstField.add(a, b)
public override fun multiply(a: MST, k: Number): MST.Binary = MstField.multiply(a, k)
public override fun scale(a: MST, value: Double): MST =
binaryOperation(GroupOperations.PLUS_OPERATION, a, number(value))
public override fun multiply(a: MST, b: MST): MST.Binary = MstField.multiply(a, b) public override fun multiply(a: MST, b: MST): MST.Binary = MstField.multiply(a, b)
public override fun divide(a: MST, b: MST): MST.Binary = MstField.divide(a, b) public override fun divide(a: MST, b: MST): MST.Binary = MstField.divide(a, b)
public override operator fun MST.unaryPlus(): MST.Unary = MstField { +this@unaryPlus } public override operator fun MST.unaryPlus(): MST.Unary = MstField { +this@unaryPlus }

@ -54,13 +54,13 @@ public inline fun <reified T : Any, A : Algebra<T>, E : Algebra<MST>> A.mst(
): MstExpression<T, A> = MstExpression(this, mstAlgebra.block()) ): MstExpression<T, A> = MstExpression(this, mstAlgebra.block())
/** /**
* Builds [MstExpression] over [Space]. * Builds [MstExpression] over [Group].
* *
* @author Alexander Nozik * @author Alexander Nozik
*/ */
public inline fun <reified T : Any, A : Space<T>> A.mstInSpace(block: MstSpace.() -> MST): MstExpression<T, A> { public inline fun <reified T : Any, A : Group<T>> A.mstInGroup(block: MstGroup.() -> MST): MstExpression<T, A> {
contract { callsInPlace(block, InvocationKind.EXACTLY_ONCE) } contract { callsInPlace(block, InvocationKind.EXACTLY_ONCE) }
return MstExpression(this, MstSpace.block()) return MstExpression(this, MstGroup.block())
} }
/** /**
@ -94,13 +94,13 @@ public inline fun <reified T : Any, A : ExtendedField<T>> A.mstInExtendedField(b
} }
/** /**
* Builds [MstExpression] over [FunctionalExpressionSpace]. * Builds [MstExpression] over [FunctionalExpressionGroup].
* *
* @author Alexander Nozik * @author Alexander Nozik
*/ */
public inline fun <reified T : Any, A : Space<T>> FunctionalExpressionSpace<T, A>.mstInSpace(block: MstSpace.() -> MST): MstExpression<T, A> { public inline fun <reified T : Any, A : Group<T>> FunctionalExpressionGroup<T, A>.mstInGroup(block: MstGroup.() -> MST): MstExpression<T, A> {
contract { callsInPlace(block, InvocationKind.EXACTLY_ONCE) } contract { callsInPlace(block, InvocationKind.EXACTLY_ONCE) }
return algebra.mstInSpace(block) return algebra.mstInGroup(block)
} }
/** /**

@ -12,12 +12,12 @@ import kotlin.test.assertEquals
internal class TestESTreeConsistencyWithInterpreter { internal class TestESTreeConsistencyWithInterpreter {
@Test @Test
fun mstSpace() { fun mstSpace() {
val res1 = MstSpace.mstInSpace { val res1 = MstGroup.mstInGroup {
binaryOperationFunction("+")( binaryOperationFunction("+")(
unaryOperationFunction("+")( unaryOperationFunction("+")(
number(3.toByte()) - (number(2.toByte()) + (multiply( number(3.toByte()) - (number(2.toByte()) + (scale(
add(number(1), number(1)), add(number(1), number(1)),
2 2.0
) + number(1.toByte()) * 3.toByte() - number(1.toByte()))) ) + number(1.toByte()) * 3.toByte() - number(1.toByte())))
), ),
@ -25,12 +25,12 @@ internal class TestESTreeConsistencyWithInterpreter {
) + bindSymbol("x") + zero ) + bindSymbol("x") + zero
}("x" to MST.Numeric(2)) }("x" to MST.Numeric(2))
val res2 = MstSpace.mstInSpace { val res2 = MstGroup.mstInGroup {
binaryOperationFunction("+")( binaryOperationFunction("+")(
unaryOperationFunction("+")( unaryOperationFunction("+")(
number(3.toByte()) - (number(2.toByte()) + (multiply( number(3.toByte()) - (number(2.toByte()) + (scale(
add(number(1), number(1)), add(number(1), number(1)),
2 2.0
) + number(1.toByte()) * 3.toByte() - number(1.toByte()))) ) + number(1.toByte()) * 3.toByte() - number(1.toByte())))
), ),
@ -46,9 +46,9 @@ internal class TestESTreeConsistencyWithInterpreter {
val res1 = ByteRing.mstInRing { val res1 = ByteRing.mstInRing {
binaryOperationFunction("+")( binaryOperationFunction("+")(
unaryOperationFunction("+")( unaryOperationFunction("+")(
(bindSymbol("x") - (2.toByte() + (multiply( (bindSymbol("x") - (2.toByte() + (scale(
add(number(1), number(1)), add(number(1), number(1)),
2 2.0
) + 1.toByte()))) * 3.0 - 1.toByte() ) + 1.toByte()))) * 3.0 - 1.toByte()
), ),
@ -59,9 +59,9 @@ internal class TestESTreeConsistencyWithInterpreter {
val res2 = ByteRing.mstInRing { val res2 = ByteRing.mstInRing {
binaryOperationFunction("+")( binaryOperationFunction("+")(
unaryOperationFunction("+")( unaryOperationFunction("+")(
(bindSymbol("x") - (2.toByte() + (multiply( (bindSymbol("x") - (2.toByte() + (scale(
add(number(1), number(1)), add(number(1), number(1)),
2 2.0
) + 1.toByte()))) * 3.0 - 1.toByte() ) + 1.toByte()))) * 3.0 - 1.toByte()
), ),
number(1) number(1)
@ -75,7 +75,7 @@ internal class TestESTreeConsistencyWithInterpreter {
fun realField() { fun realField() {
val res1 = RealField.mstInField { val res1 = RealField.mstInField {
+(3 - 2 + 2 * number(1) + 1.0) + binaryOperationFunction("+")( +(3 - 2 + 2 * number(1) + 1.0) + binaryOperationFunction("+")(
(3.0 - (bindSymbol("x") + (multiply(add(number(1.0), number(1.0)), 2) + 1.0))) * 3 - 1.0 (3.0 - (bindSymbol("x") + (scale(add(number(1.0), number(1.0)), 2.0) + 1.0))) * 3 - 1.0
+ number(1), + number(1),
number(1) / 2 + number(2.0) * one number(1) / 2 + number(2.0) * one
) + zero ) + zero
@ -83,7 +83,7 @@ internal class TestESTreeConsistencyWithInterpreter {
val res2 = RealField.mstInField { val res2 = RealField.mstInField {
+(3 - 2 + 2 * number(1) + 1.0) + binaryOperationFunction("+")( +(3 - 2 + 2 * number(1) + 1.0) + binaryOperationFunction("+")(
(3.0 - (bindSymbol("x") + (multiply(add(number(1.0), number(1.0)), 2) + 1.0))) * 3 - 1.0 (3.0 - (bindSymbol("x") + (scale(add(number(1.0), number(1.0)), 2.0) + 1.0))) * 3 - 1.0
+ number(1), + number(1),
number(1) / 2 + number(2.0) * one number(1) / 2 + number(2.0) * one
) + zero ) + zero
@ -96,7 +96,7 @@ internal class TestESTreeConsistencyWithInterpreter {
fun complexField() { fun complexField() {
val res1 = ComplexField.mstInField { val res1 = ComplexField.mstInField {
+(3 - 2 + 2 * number(1) + 1.0) + binaryOperationFunction("+")( +(3 - 2 + 2 * number(1) + 1.0) + binaryOperationFunction("+")(
(3.0 - (bindSymbol("x") + (multiply(add(number(1.0), number(1.0)), 2) + 1.0))) * 3 - 1.0 (3.0 - (bindSymbol("x") + (scale(add(number(1.0), number(1.0)), 2.0) + 1.0))) * 3 - 1.0
+ number(1), + number(1),
number(1) / 2 + number(2.0) * one number(1) / 2 + number(2.0) * one
) + zero ) + zero
@ -104,7 +104,7 @@ internal class TestESTreeConsistencyWithInterpreter {
val res2 = ComplexField.mstInField { val res2 = ComplexField.mstInField {
+(3 - 2 + 2 * number(1) + 1.0) + binaryOperationFunction("+")( +(3 - 2 + 2 * number(1) + 1.0) + binaryOperationFunction("+")(
(3.0 - (bindSymbol("x") + (multiply(add(number(1.0), number(1.0)), 2) + 1.0))) * 3 - 1.0 (3.0 - (bindSymbol("x") + (scale(add(number(1.0), number(1.0)), 2.0) + 1.0))) * 3 - 1.0
+ number(1), + number(1),
number(1) / 2 + number(2.0) * one number(1) / 2 + number(2.0) * one
) + zero ) + zero

@ -2,7 +2,7 @@ package space.kscience.kmath.estree
import space.kscience.kmath.ast.mstInExtendedField import space.kscience.kmath.ast.mstInExtendedField
import space.kscience.kmath.ast.mstInField import space.kscience.kmath.ast.mstInField
import space.kscience.kmath.ast.mstInSpace import space.kscience.kmath.ast.mstInGroup
import space.kscience.kmath.expressions.invoke import space.kscience.kmath.expressions.invoke
import space.kscience.kmath.operations.RealField import space.kscience.kmath.operations.RealField
import kotlin.random.Random import kotlin.random.Random
@ -12,14 +12,14 @@ import kotlin.test.assertEquals
internal class TestESTreeOperationsSupport { internal class TestESTreeOperationsSupport {
@Test @Test
fun testUnaryOperationInvocation() { fun testUnaryOperationInvocation() {
val expression = RealField.mstInSpace { -bindSymbol("x") }.compile() val expression = RealField.mstInGroup { -bindSymbol("x") }.compile()
val res = expression("x" to 2.0) val res = expression("x" to 2.0)
assertEquals(-2.0, res) assertEquals(-2.0, res)
} }
@Test @Test
fun testBinaryOperationInvocation() { fun testBinaryOperationInvocation() {
val expression = RealField.mstInSpace { -bindSymbol("x") + number(1.0) }.compile() val expression = RealField.mstInGroup { -bindSymbol("x") + number(1.0) }.compile()
val res = expression("x" to 2.0) val res = expression("x" to 2.0)
assertEquals(-1.0, res) assertEquals(-1.0, res)
} }

@ -14,9 +14,9 @@ import com.github.h0tk3y.betterParse.lexer.regexToken
import com.github.h0tk3y.betterParse.parser.ParseResult import com.github.h0tk3y.betterParse.parser.ParseResult
import com.github.h0tk3y.betterParse.parser.Parser import com.github.h0tk3y.betterParse.parser.Parser
import space.kscience.kmath.operations.FieldOperations import space.kscience.kmath.operations.FieldOperations
import space.kscience.kmath.operations.GroupOperations
import space.kscience.kmath.operations.PowerOperations import space.kscience.kmath.operations.PowerOperations
import space.kscience.kmath.operations.RingOperations import space.kscience.kmath.operations.RingOperations
import space.kscience.kmath.operations.SpaceOperations
/** /**
* better-parse implementation of grammar defined in the ArithmeticsEvaluator.g4. * better-parse implementation of grammar defined in the ArithmeticsEvaluator.g4.
@ -55,7 +55,7 @@ public object ArithmeticsEvaluator : Grammar<MST>() {
.or(binaryFunction) .or(binaryFunction)
.or(unaryFunction) .or(unaryFunction)
.or(singular) .or(singular)
.or(-minus and parser(ArithmeticsEvaluator::term) map { MST.Unary(SpaceOperations.MINUS_OPERATION, it) }) .or(-minus and parser(ArithmeticsEvaluator::term) map { MST.Unary(GroupOperations.MINUS_OPERATION, it) })
.or(-lpar and parser(ArithmeticsEvaluator::subSumChain) and -rpar) .or(-lpar and parser(ArithmeticsEvaluator::subSumChain) and -rpar)
private val powChain: Parser<MST> by leftAssociative(term = term, operator = pow) { a, _, b -> private val powChain: Parser<MST> by leftAssociative(term = term, operator = pow) { a, _, b ->
@ -77,9 +77,9 @@ public object ArithmeticsEvaluator : Grammar<MST>() {
operator = plus or minus use TokenMatch::type operator = plus or minus use TokenMatch::type
) { a, op, b -> ) { a, op, b ->
if (op == plus) if (op == plus)
MST.Binary(SpaceOperations.PLUS_OPERATION, a, b) MST.Binary(GroupOperations.PLUS_OPERATION, a, b)
else else
MST.Binary(SpaceOperations.MINUS_OPERATION, a, b) MST.Binary(GroupOperations.MINUS_OPERATION, a, b)
} }
override val rootParser: Parser<MST> by subSumChain override val rootParser: Parser<MST> by subSumChain

@ -12,12 +12,12 @@ import kotlin.test.assertEquals
internal class TestAsmConsistencyWithInterpreter { internal class TestAsmConsistencyWithInterpreter {
@Test @Test
fun mstSpace() { fun mstSpace() {
val res1 = MstSpace.mstInSpace { val res1 = MstGroup.mstInGroup {
binaryOperationFunction("+")( binaryOperationFunction("+")(
unaryOperationFunction("+")( unaryOperationFunction("+")(
number(3.toByte()) - (number(2.toByte()) + (multiply( number(3.toByte()) - (number(2.toByte()) + (scale(
add(number(1), number(1)), add(number(1), number(1)),
2 2.0
) + number(1.toByte()) * 3.toByte() - number(1.toByte()))) ) + number(1.toByte()) * 3.toByte() - number(1.toByte())))
), ),
@ -25,12 +25,12 @@ internal class TestAsmConsistencyWithInterpreter {
) + bindSymbol("x") + zero ) + bindSymbol("x") + zero
}("x" to MST.Numeric(2)) }("x" to MST.Numeric(2))
val res2 = MstSpace.mstInSpace { val res2 = MstGroup.mstInGroup {
binaryOperationFunction("+")( binaryOperationFunction("+")(
unaryOperationFunction("+")( unaryOperationFunction("+")(
number(3.toByte()) - (number(2.toByte()) + (multiply( number(3.toByte()) - (number(2.toByte()) + (scale(
add(number(1), number(1)), add(number(1), number(1)),
2 2.0
) + number(1.toByte()) * 3.toByte() - number(1.toByte()))) ) + number(1.toByte()) * 3.toByte() - number(1.toByte())))
), ),
@ -46,9 +46,9 @@ internal class TestAsmConsistencyWithInterpreter {
val res1 = ByteRing.mstInRing { val res1 = ByteRing.mstInRing {
binaryOperationFunction("+")( binaryOperationFunction("+")(
unaryOperationFunction("+")( unaryOperationFunction("+")(
(bindSymbol("x") - (2.toByte() + (multiply( (bindSymbol("x") - (2.toByte() + (scale(
add(number(1), number(1)), add(number(1), number(1)),
2 2.0
) + 1.toByte()))) * 3.0 - 1.toByte() ) + 1.toByte()))) * 3.0 - 1.toByte()
), ),
@ -59,9 +59,9 @@ internal class TestAsmConsistencyWithInterpreter {
val res2 = ByteRing.mstInRing { val res2 = ByteRing.mstInRing {
binaryOperationFunction("+")( binaryOperationFunction("+")(
unaryOperationFunction("+")( unaryOperationFunction("+")(
(bindSymbol("x") - (2.toByte() + (multiply( (bindSymbol("x") - (2.toByte() + (scale(
add(number(1), number(1)), add(number(1), number(1)),
2 2.0
) + 1.toByte()))) * 3.0 - 1.toByte() ) + 1.toByte()))) * 3.0 - 1.toByte()
), ),
number(1) number(1)
@ -75,7 +75,7 @@ internal class TestAsmConsistencyWithInterpreter {
fun realField() { fun realField() {
val res1 = RealField.mstInField { val res1 = RealField.mstInField {
+(3 - 2 + 2 * number(1) + 1.0) + binaryOperationFunction("+")( +(3 - 2 + 2 * number(1) + 1.0) + binaryOperationFunction("+")(
(3.0 - (bindSymbol("x") + (multiply(add(number(1.0), number(1.0)), 2) + 1.0))) * 3 - 1.0 (3.0 - (bindSymbol("x") + (scale(add(number(1.0), number(1.0)), 2.0) + 1.0))) * 3 - 1.0
+ number(1), + number(1),
number(1) / 2 + number(2.0) * one number(1) / 2 + number(2.0) * one
) + zero ) + zero
@ -83,7 +83,7 @@ internal class TestAsmConsistencyWithInterpreter {
val res2 = RealField.mstInField { val res2 = RealField.mstInField {
+(3 - 2 + 2 * number(1) + 1.0) + binaryOperationFunction("+")( +(3 - 2 + 2 * number(1) + 1.0) + binaryOperationFunction("+")(
(3.0 - (bindSymbol("x") + (multiply(add(number(1.0), number(1.0)), 2) + 1.0))) * 3 - 1.0 (3.0 - (bindSymbol("x") + (scale(add(number(1.0), number(1.0)), 2.0) + 1.0))) * 3 - 1.0
+ number(1), + number(1),
number(1) / 2 + number(2.0) * one number(1) / 2 + number(2.0) * one
) + zero ) + zero
@ -96,7 +96,7 @@ internal class TestAsmConsistencyWithInterpreter {
fun complexField() { fun complexField() {
val res1 = ComplexField.mstInField { val res1 = ComplexField.mstInField {
+(3 - 2 + 2 * number(1) + 1.0) + binaryOperationFunction("+")( +(3 - 2 + 2 * number(1) + 1.0) + binaryOperationFunction("+")(
(3.0 - (bindSymbol("x") + (multiply(add(number(1.0), number(1.0)), 2) + 1.0))) * 3 - 1.0 (3.0 - (bindSymbol("x") + (scale(add(number(1.0), number(1.0)), 2.0) + 1.0))) * 3 - 1.0
+ number(1), + number(1),
number(1) / 2 + number(2.0) * one number(1) / 2 + number(2.0) * one
) + zero ) + zero
@ -104,7 +104,7 @@ internal class TestAsmConsistencyWithInterpreter {
val res2 = ComplexField.mstInField { val res2 = ComplexField.mstInField {
+(3 - 2 + 2 * number(1) + 1.0) + binaryOperationFunction("+")( +(3 - 2 + 2 * number(1) + 1.0) + binaryOperationFunction("+")(
(3.0 - (bindSymbol("x") + (multiply(add(number(1.0), number(1.0)), 2) + 1.0))) * 3 - 1.0 (3.0 - (bindSymbol("x") + (scale(add(number(1.0), number(1.0)), 2.0) + 1.0))) * 3 - 1.0
+ number(1), + number(1),
number(1) / 2 + number(2.0) * one number(1) / 2 + number(2.0) * one
) + zero ) + zero

@ -2,7 +2,7 @@ package space.kscience.kmath.asm
import space.kscience.kmath.ast.mstInExtendedField import space.kscience.kmath.ast.mstInExtendedField
import space.kscience.kmath.ast.mstInField import space.kscience.kmath.ast.mstInField
import space.kscience.kmath.ast.mstInSpace import space.kscience.kmath.ast.mstInGroup
import space.kscience.kmath.expressions.invoke import space.kscience.kmath.expressions.invoke
import space.kscience.kmath.operations.RealField import space.kscience.kmath.operations.RealField
import kotlin.random.Random import kotlin.random.Random
@ -12,14 +12,14 @@ import kotlin.test.assertEquals
internal class TestAsmOperationsSupport { internal class TestAsmOperationsSupport {
@Test @Test
fun testUnaryOperationInvocation() { fun testUnaryOperationInvocation() {
val expression = RealField.mstInSpace { -bindSymbol("x") }.compile() val expression = RealField.mstInGroup { -bindSymbol("x") }.compile()
val res = expression("x" to 2.0) val res = expression("x" to 2.0)
assertEquals(-2.0, res) assertEquals(-2.0, res)
} }
@Test @Test
fun testBinaryOperationInvocation() { fun testBinaryOperationInvocation() {
val expression = RealField.mstInSpace { -bindSymbol("x") + number(1.0) }.compile() val expression = RealField.mstInGroup { -bindSymbol("x") + number(1.0) }.compile()
val res = expression("x" to 2.0) val res = expression("x" to 2.0)
assertEquals(-1.0, res) assertEquals(-1.0, res)
} }

@ -4,7 +4,7 @@ import org.apache.commons.math3.analysis.differentiation.DerivativeStructure
import space.kscience.kmath.expressions.* import space.kscience.kmath.expressions.*
import space.kscience.kmath.misc.UnstableKMathAPI import space.kscience.kmath.misc.UnstableKMathAPI
import space.kscience.kmath.operations.ExtendedField import space.kscience.kmath.operations.ExtendedField
import space.kscience.kmath.operations.RingWithNumbers import space.kscience.kmath.operations.NumbersAddOperations
/** /**
* A field over commons-math [DerivativeStructure]. * A field over commons-math [DerivativeStructure].
@ -16,7 +16,8 @@ import space.kscience.kmath.operations.RingWithNumbers
public class DerivativeStructureField( public class DerivativeStructureField(
public val order: Int, public val order: Int,
bindings: Map<Symbol, Double>, bindings: Map<Symbol, Double>,
) : ExtendedField<DerivativeStructure>, ExpressionAlgebra<Double, DerivativeStructure>, RingWithNumbers<DerivativeStructure> { ) : ExtendedField<DerivativeStructure>, ExpressionAlgebra<Double, DerivativeStructure>,
NumbersAddOperations<DerivativeStructure> {
public val numberOfVariables: Int = bindings.size public val numberOfVariables: Int = bindings.size
public override val zero: DerivativeStructure by lazy { DerivativeStructure(numberOfVariables, order) } public override val zero: DerivativeStructure by lazy { DerivativeStructure(numberOfVariables, order) }
@ -62,13 +63,11 @@ public class DerivativeStructureField(
public fun DerivativeStructure.derivative(vararg symbols: Symbol): Double = derivative(symbols.toList()) public fun DerivativeStructure.derivative(vararg symbols: Symbol): Double = derivative(symbols.toList())
override fun DerivativeStructure.unaryMinus(): DerivativeStructure = negate()
public override fun add(a: DerivativeStructure, b: DerivativeStructure): DerivativeStructure = a.add(b) public override fun add(a: DerivativeStructure, b: DerivativeStructure): DerivativeStructure = a.add(b)
public override fun multiply(a: DerivativeStructure, k: Number): DerivativeStructure = when (k) { public override fun scale(a: DerivativeStructure, value: Double): DerivativeStructure = a.multiply(value)
is Double -> a.multiply(k)
is Int -> a.multiply(k)
else -> a.multiply(k.toDouble())
}
public override fun multiply(a: DerivativeStructure, b: DerivativeStructure): DerivativeStructure = a.multiply(b) public override fun multiply(a: DerivativeStructure, b: DerivativeStructure): DerivativeStructure = a.multiply(b)
public override fun divide(a: DerivativeStructure, b: DerivativeStructure): DerivativeStructure = a.divide(b) public override fun divide(a: DerivativeStructure, b: DerivativeStructure): DerivativeStructure = a.divide(b)

@ -3,11 +3,13 @@ package space.kscience.kmath.commons.linear
import org.apache.commons.math3.linear.* import org.apache.commons.math3.linear.*
import space.kscience.kmath.linear.* import space.kscience.kmath.linear.*
import space.kscience.kmath.misc.UnstableKMathAPI import space.kscience.kmath.misc.UnstableKMathAPI
import space.kscience.kmath.nd.NDStructure
import space.kscience.kmath.operations.RealField
import space.kscience.kmath.structures.RealBuffer import space.kscience.kmath.structures.RealBuffer
import kotlin.reflect.KClass import kotlin.reflect.KClass
import kotlin.reflect.cast import kotlin.reflect.cast
public inline class CMMatrix(public val origin: RealMatrix) : Matrix<Double> { public class CMMatrix(public val origin: RealMatrix) : Matrix<Double> {
public override val rowNum: Int get() = origin.rowDimension public override val rowNum: Int get() = origin.rowDimension
public override val colNum: Int get() = origin.columnDimension public override val colNum: Int get() = origin.columnDimension
@ -50,12 +52,17 @@ public inline class CMMatrix(public val origin: RealMatrix) : Matrix<Double> {
}?.let(type::cast) }?.let(type::cast)
public override operator fun get(i: Int, j: Int): Double = origin.getEntry(i, j) public override operator fun get(i: Int, j: Int): Double = origin.getEntry(i, j)
override fun equals(other: Any?): Boolean {
if (this === other) return true
if (other !is NDStructure<*>) return false
return NDStructure.contentEquals(this, other)
}
override fun hashCode(): Int = origin.hashCode()
} }
public inline class CMVector(public val origin: RealVector) : Point<Double> {
public fun RealMatrix.asMatrix(): CMMatrix = CMMatrix(this)
public class CMVector(public val origin: RealVector) : Point<Double> {
public override val size: Int get() = origin.dimension public override val size: Int get() = origin.dimension
public override operator fun get(index: Int): Double = origin.getEntry(index) public override operator fun get(index: Int): Double = origin.getEntry(index)
@ -63,16 +70,17 @@ public class CMVector(public val origin: RealVector) : Point<Double> {
public override operator fun iterator(): Iterator<Double> = origin.toArray().iterator() public override operator fun iterator(): Iterator<Double> = origin.toArray().iterator()
} }
public fun Point<Double>.toCM(): CMVector = if (this is CMVector) this else {
val array = DoubleArray(size) { this[it] }
CMVector(ArrayRealVector(array))
}
public fun RealVector.toPoint(): CMVector = CMVector(this) public fun RealVector.toPoint(): CMVector = CMVector(this)
public object CMMatrixContext : MatrixContext<Double, CMMatrix> { public object CMLinearSpace : LinearSpace<Double, RealField> {
public override fun produce(rows: Int, columns: Int, initializer: (i: Int, j: Int) -> Double): CMMatrix { override val elementAlgebra: RealField get() = RealField
val array = Array(rows) { i -> DoubleArray(columns) { j -> initializer(i, j) } }
public override fun buildMatrix(
rows: Int,
columns: Int,
initializer: RealField.(i: Int, j: Int) -> Double,
): CMMatrix {
val array = Array(rows) { i -> DoubleArray(columns) { j -> RealField.initializer(i, j) } }
return CMMatrix(Array2DRowRealMatrix(array)) return CMMatrix(Array2DRowRealMatrix(array))
} }
@ -82,30 +90,50 @@ public object CMMatrixContext : MatrixContext<Double, CMMatrix> {
else -> { else -> {
//TODO add feature analysis //TODO add feature analysis
val array = Array(rowNum) { i -> DoubleArray(colNum) { j -> get(i, j) } } val array = Array(rowNum) { i -> DoubleArray(colNum) { j -> get(i, j) } }
CMMatrix(Array2DRowRealMatrix(array)) Array2DRowRealMatrix(array).wrap()
} }
} }
public fun Point<Double>.toCM(): CMVector = if (this is CMVector) this else {
val array = DoubleArray(size) { this[it] }
ArrayRealVector(array).wrap()
}
internal fun RealMatrix.wrap(): CMMatrix = CMMatrix(this)
internal fun RealVector.wrap(): CMVector = CMVector(this)
override fun buildVector(size: Int, initializer: RealField.(Int) -> Double): Point<Double> =
ArrayRealVector(DoubleArray(size) { RealField.initializer(it) }).wrap()
override fun Matrix<Double>.plus(other: Matrix<Double>): CMMatrix =
toCM().origin.add(other.toCM().origin).wrap()
override fun Point<Double>.plus(other: Point<Double>): CMVector =
toCM().origin.add(other.toCM().origin).wrap()
override fun Point<Double>.minus(other: Point<Double>): CMVector =
toCM().origin.subtract(other.toCM().origin).wrap()
public override fun Matrix<Double>.dot(other: Matrix<Double>): CMMatrix = public override fun Matrix<Double>.dot(other: Matrix<Double>): CMMatrix =
CMMatrix(toCM().origin.multiply(other.toCM().origin)) toCM().origin.multiply(other.toCM().origin).wrap()
public override fun Matrix<Double>.dot(vector: Point<Double>): CMVector = public override fun Matrix<Double>.dot(vector: Point<Double>): CMVector =
CMVector(toCM().origin.preMultiply(vector.toCM().origin)) toCM().origin.preMultiply(vector.toCM().origin).wrap()
public override operator fun Matrix<Double>.unaryMinus(): CMMatrix = public override operator fun Matrix<Double>.minus(other: Matrix<Double>): CMMatrix =
produce(rowNum, colNum) { i, j -> -get(i, j) } toCM().origin.subtract(other.toCM().origin).wrap()
public override fun add(a: Matrix<Double>, b: Matrix<Double>): CMMatrix =
CMMatrix(a.toCM().origin.multiply(b.toCM().origin))
public override operator fun Matrix<Double>.minus(b: Matrix<Double>): CMMatrix =
CMMatrix(toCM().origin.subtract(b.toCM().origin))
public override fun multiply(a: Matrix<Double>, k: Number): CMMatrix =
CMMatrix(a.toCM().origin.scalarMultiply(k.toDouble()))
public override operator fun Matrix<Double>.times(value: Double): CMMatrix = public override operator fun Matrix<Double>.times(value: Double): CMMatrix =
produce(rowNum, colNum) { i, j -> get(i, j) * value } toCM().origin.scalarMultiply(value).wrap()
override fun Double.times(m: Matrix<Double>): CMMatrix =
m * this
override fun Point<Double>.times(value: Double): CMVector =
toCM().origin.mapMultiply(value).wrap()
override fun Double.times(v: Point<Double>): CMVector =
v * this
} }
public operator fun CMMatrix.plus(other: CMMatrix): CMMatrix = public operator fun CMMatrix.plus(other: CMMatrix): CMMatrix =

@ -12,7 +12,7 @@ public enum class CMDecomposition {
CHOLESKY CHOLESKY
} }
public fun CMMatrixContext.solver( public fun CMLinearSpace.solver(
a: Matrix<Double>, a: Matrix<Double>,
decomposition: CMDecomposition = CMDecomposition.LUP decomposition: CMDecomposition = CMDecomposition.LUP
): DecompositionSolver = when (decomposition) { ): DecompositionSolver = when (decomposition) {
@ -23,19 +23,19 @@ public fun CMMatrixContext.solver(
CMDecomposition.CHOLESKY -> CholeskyDecomposition(a.toCM().origin).solver CMDecomposition.CHOLESKY -> CholeskyDecomposition(a.toCM().origin).solver
} }
public fun CMMatrixContext.solve( public fun CMLinearSpace.solve(
a: Matrix<Double>, a: Matrix<Double>,
b: Matrix<Double>, b: Matrix<Double>,
decomposition: CMDecomposition = CMDecomposition.LUP decomposition: CMDecomposition = CMDecomposition.LUP
): CMMatrix = solver(a, decomposition).solve(b.toCM().origin).asMatrix() ): CMMatrix = solver(a, decomposition).solve(b.toCM().origin).wrap()
public fun CMMatrixContext.solve( public fun CMLinearSpace.solve(
a: Matrix<Double>, a: Matrix<Double>,
b: Point<Double>, b: Point<Double>,
decomposition: CMDecomposition = CMDecomposition.LUP decomposition: CMDecomposition = CMDecomposition.LUP
): CMVector = solver(a, decomposition).solve(b.toCM().origin).toPoint() ): CMVector = solver(a, decomposition).solve(b.toCM().origin).toPoint()
public fun CMMatrixContext.inverse( public fun CMLinearSpace.inverse(
a: Matrix<Double>, a: Matrix<Double>,
decomposition: CMDecomposition = CMDecomposition.LUP decomposition: CMDecomposition = CMDecomposition.LUP
): CMMatrix = solver(a, decomposition).inverse.asMatrix() ): CMMatrix = solver(a, decomposition).inverse.wrap()

@ -4,10 +4,7 @@ import space.kscience.kmath.memory.MemoryReader
import space.kscience.kmath.memory.MemorySpec import space.kscience.kmath.memory.MemorySpec
import space.kscience.kmath.memory.MemoryWriter import space.kscience.kmath.memory.MemoryWriter
import space.kscience.kmath.misc.UnstableKMathAPI import space.kscience.kmath.misc.UnstableKMathAPI
import space.kscience.kmath.operations.ExtendedField import space.kscience.kmath.operations.*
import space.kscience.kmath.operations.FieldElement
import space.kscience.kmath.operations.Norm
import space.kscience.kmath.operations.RingWithNumbers
import space.kscience.kmath.structures.Buffer import space.kscience.kmath.structures.Buffer
import space.kscience.kmath.structures.MemoryBuffer import space.kscience.kmath.structures.MemoryBuffer
import space.kscience.kmath.structures.MutableBuffer import space.kscience.kmath.structures.MutableBuffer
@ -47,7 +44,8 @@ private val PI_DIV_2 = Complex(PI / 2, 0)
* A field of [Complex]. * A field of [Complex].
*/ */
@OptIn(UnstableKMathAPI::class) @OptIn(UnstableKMathAPI::class)
public object ComplexField : ExtendedField<Complex>, Norm<Complex, Complex>, RingWithNumbers<Complex> { public object ComplexField : ExtendedField<Complex>, Norm<Complex, Complex>, NumbersAddOperations<Complex>,
ScaleOperations<Complex> {
public override val zero: Complex = 0.0.toComplex() public override val zero: Complex = 0.0.toComplex()
public override val one: Complex = 1.0.toComplex() public override val one: Complex = 1.0.toComplex()
@ -56,8 +54,14 @@ public object ComplexField : ExtendedField<Complex>, Norm<Complex, Complex>, Rin
*/ */
public val i: Complex by lazy { Complex(0.0, 1.0) } public val i: Complex by lazy { Complex(0.0, 1.0) }
override fun Complex.unaryMinus(): Complex = Complex(-re, -im)
override fun number(value: Number): Complex = Complex(value.toDouble(), 0.0)
override fun scale(a: Complex, value: Double): Complex = Complex(a.re * value, a.im * value)
public override fun add(a: Complex, b: Complex): Complex = Complex(a.re + b.re, a.im + b.im) public override fun add(a: Complex, b: Complex): Complex = Complex(a.re + b.re, a.im + b.im)
public override fun multiply(a: Complex, k: Number): Complex = Complex(a.re * k.toDouble(), a.im * k.toDouble()) // public override fun multiply(a: Complex, k: Number): Complex = Complex(a.re * k.toDouble(), a.im * k.toDouble())
public override fun multiply(a: Complex, b: Complex): Complex = public override fun multiply(a: Complex, b: Complex): Complex =
Complex(a.re * b.re - a.im * b.im, a.re * b.im + a.im * b.re) Complex(a.re * b.re - a.im * b.im, a.re * b.im + a.im * b.re)
@ -86,8 +90,10 @@ public object ComplexField : ExtendedField<Complex>, Norm<Complex, Complex>, Rin
} }
} }
public override fun sin(arg: Complex): Complex = i * (exp(-i * arg) - exp(i * arg)) / 2 override operator fun Complex.div(k: Number): Complex = Complex(re / k.toDouble(), im / k.toDouble())
public override fun cos(arg: Complex): Complex = (exp(-i * arg) + exp(i * arg)) / 2
public override fun sin(arg: Complex): Complex = i * (exp(-i * arg) - exp(i * arg)) / 2.0
public override fun cos(arg: Complex): Complex = (exp(-i * arg) + exp(i * arg)) / 2.0
public override fun tan(arg: Complex): Complex { public override fun tan(arg: Complex): Complex {
val e1 = exp(-i * arg) val e1 = exp(-i * arg)
@ -159,7 +165,8 @@ public object ComplexField : ExtendedField<Complex>, Norm<Complex, Complex>, Rin
public override fun norm(arg: Complex): Complex = sqrt(arg.conjugate * arg) public override fun norm(arg: Complex): Complex = sqrt(arg.conjugate * arg)
public override fun bindSymbol(value: String): Complex = if (value == "i") i else super<ExtendedField>.bindSymbol(value) public override fun bindSymbol(value: String): Complex =
if (value == "i") i else super<ExtendedField>.bindSymbol(value)
} }
/** /**
@ -181,7 +188,8 @@ public data class Complex(val re: Double, val im: Double) : FieldElement<Complex
public override val objectSize: Int public override val objectSize: Int
get() = 16 get() = 16
public override fun MemoryReader.read(offset: Int): Complex = Complex(readDouble(offset), readDouble(offset + 8)) public override fun MemoryReader.read(offset: Int): Complex =
Complex(readDouble(offset), readDouble(offset + 8))
public override fun MemoryWriter.write(offset: Int, value: Complex) { public override fun MemoryWriter.write(offset: Int, value: Complex) {
writeDouble(offset, value.re) writeDouble(offset, value.re)

@ -6,7 +6,7 @@ import space.kscience.kmath.nd.NDAlgebra
import space.kscience.kmath.nd.NDBuffer import space.kscience.kmath.nd.NDBuffer
import space.kscience.kmath.nd.NDStructure import space.kscience.kmath.nd.NDStructure
import space.kscience.kmath.operations.ExtendedField import space.kscience.kmath.operations.ExtendedField
import space.kscience.kmath.operations.RingWithNumbers import space.kscience.kmath.operations.NumbersAddOperations
import space.kscience.kmath.structures.Buffer import space.kscience.kmath.structures.Buffer
import kotlin.contracts.InvocationKind import kotlin.contracts.InvocationKind
import kotlin.contracts.contract import kotlin.contracts.contract
@ -19,7 +19,7 @@ import kotlin.contracts.contract
public class ComplexNDField( public class ComplexNDField(
shape: IntArray, shape: IntArray,
) : BufferedNDField<Complex, ComplexField>(shape, ComplexField, Buffer.Companion::complex), ) : BufferedNDField<Complex, ComplexField>(shape, ComplexField, Buffer.Companion::complex),
RingWithNumbers<NDStructure<Complex>>, NumbersAddOperations<NDStructure<Complex>>,
ExtendedField<NDStructure<Complex>> { ExtendedField<NDStructure<Complex>> {
override val zero: NDBuffer<Complex> by lazy { produce { zero } } override val zero: NDBuffer<Complex> by lazy { produce { zero } }
@ -29,6 +29,7 @@ public class ComplexNDField(
val d = value.toComplex() // minimize conversions val d = value.toComplex() // minimize conversions
return produce { d } return produce { d }
} }
// //
// @Suppress("OVERRIDE_BY_INLINE") // @Suppress("OVERRIDE_BY_INLINE")
// override inline fun map( // override inline fun map(
@ -75,25 +76,25 @@ public class ComplexNDField(
// return BufferedNDFieldElement(this, buffer) // return BufferedNDFieldElement(this, buffer)
// } // }
override fun power(arg: NDStructure<Complex>, pow: Number): NDBuffer<Complex> = arg.map() { power(it, pow) } override fun power(arg: NDStructure<Complex>, pow: Number): NDBuffer<Complex> = arg.map { power(it, pow) }
override fun exp(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map() { exp(it) } override fun exp(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map { exp(it) }
override fun ln(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map() { ln(it) } override fun ln(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map { ln(it) }
override fun sin(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map() { sin(it) } override fun sin(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map { sin(it) }
override fun cos(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map() { cos(it) } override fun cos(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map { cos(it) }
override fun tan(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map() { tan(it) } override fun tan(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map { tan(it) }
override fun asin(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map() { asin(it) } override fun asin(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map { asin(it) }
override fun acos(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map() { acos(it) } override fun acos(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map { acos(it) }
override fun atan(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map() { atan(it) } override fun atan(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map { atan(it) }
override fun sinh(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map() { sinh(it) } override fun sinh(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map { sinh(it) }
override fun cosh(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map() { cosh(it) } override fun cosh(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map { cosh(it) }
override fun tanh(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map() { tanh(it) } override fun tanh(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map { tanh(it) }
override fun asinh(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map() { asinh(it) } override fun asinh(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map { asinh(it) }
override fun acosh(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map() { acosh(it) } override fun acosh(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map { acosh(it) }
override fun atanh(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map() { atanh(it) } override fun atanh(arg: NDStructure<Complex>): NDBuffer<Complex> = arg.map { atanh(it) }
} }

@ -37,7 +37,7 @@ public val Quaternion.r: Double
*/ */
@OptIn(UnstableKMathAPI::class) @OptIn(UnstableKMathAPI::class)
public object QuaternionField : Field<Quaternion>, Norm<Quaternion, Quaternion>, PowerOperations<Quaternion>, public object QuaternionField : Field<Quaternion>, Norm<Quaternion, Quaternion>, PowerOperations<Quaternion>,
ExponentialOperations<Quaternion>, RingWithNumbers<Quaternion> { ExponentialOperations<Quaternion>, NumbersAddOperations<Quaternion>, ScaleOperations<Quaternion> {
override val zero: Quaternion = 0.toQuaternion() override val zero: Quaternion = 0.toQuaternion()
override val one: Quaternion = 1.toQuaternion() override val one: Quaternion = 1.toQuaternion()
@ -59,10 +59,8 @@ public object QuaternionField : Field<Quaternion>, Norm<Quaternion, Quaternion>,
public override fun add(a: Quaternion, b: Quaternion): Quaternion = public override fun add(a: Quaternion, b: Quaternion): Quaternion =
Quaternion(a.w + b.w, a.x + b.x, a.y + b.y, a.z + b.z) Quaternion(a.w + b.w, a.x + b.x, a.y + b.y, a.z + b.z)
public override fun multiply(a: Quaternion, k: Number): Quaternion { public override fun scale(a: Quaternion, value: Double): Quaternion =
val d = k.toDouble() Quaternion(a.w * value, a.x * value, a.y * value, a.z * value)
return Quaternion(a.w * d, a.x * d, a.y * d, a.z * d)
}
public override fun multiply(a: Quaternion, b: Quaternion): Quaternion = Quaternion( public override fun multiply(a: Quaternion, b: Quaternion): Quaternion = Quaternion(
a.w * b.w - a.x * b.x - a.y * b.y - a.z * b.z, a.w * b.w - a.x * b.x - a.y * b.y - a.z * b.z,
@ -173,6 +171,15 @@ public object QuaternionField : Field<Quaternion>, Norm<Quaternion, Quaternion>,
"k" -> k "k" -> k
else -> super<Field>.bindSymbol(value) else -> super<Field>.bindSymbol(value)
} }
override fun number(value: Number): Quaternion =value.toQuaternion()
public override fun sinh(arg: Quaternion): Quaternion = (exp(arg) - exp(-arg)) / 2.0
public override fun cosh(arg: Quaternion): Quaternion = (exp(arg) + exp(-arg)) / 2.0
public override fun tanh(arg: Quaternion): Quaternion = (exp(arg) - exp(-arg)) / (exp(-arg) + exp(arg))
public override fun asinh(arg: Quaternion): Quaternion = ln(sqrt(arg * arg + one) + arg)
public override fun acosh(arg: Quaternion): Quaternion = ln(arg + sqrt((arg - one) * (arg + one)))
public override fun atanh(arg: Quaternion): Quaternion = (ln(arg + one) - ln(one - arg)) / 2.0
} }
/** /**
@ -207,8 +214,7 @@ public data class Quaternion(
require(!z.isNaN()) { "x-component of quaternion is not-a-number" } require(!z.isNaN()) { "x-component of quaternion is not-a-number" }
} }
public override val context: QuaternionField public override val context: QuaternionField get() = QuaternionField
get() = QuaternionField
/** /**
* Returns a string representation of this quaternion. * Returns a string representation of this quaternion.

@ -4,7 +4,6 @@ import space.kscience.kmath.expressions.FunctionalExpressionField
import space.kscience.kmath.expressions.bindSymbol import space.kscience.kmath.expressions.bindSymbol
import space.kscience.kmath.expressions.invoke import space.kscience.kmath.expressions.invoke
import space.kscience.kmath.expressions.symbol import space.kscience.kmath.expressions.symbol
import space.kscience.kmath.operations.invoke
import kotlin.test.Test import kotlin.test.Test
import kotlin.test.assertEquals import kotlin.test.assertEquals
@ -13,14 +12,11 @@ internal class ExpressionFieldForComplexTest {
@Test @Test
fun testComplex() { fun testComplex() {
val context = FunctionalExpressionField(ComplexField) val expression = FunctionalExpressionField(ComplexField).run {
val expression = context {
val x = bindSymbol(x) val x = bindSymbol(x)
x * x + 2 * x + one x * x + 2 * x + one
} }
assertEquals(expression(x to Complex(1.0, 0.0)), Complex(4.0, 0.0)) assertEquals(expression(x to Complex(1.0, 0.0)), Complex(4.0, 0.0))
//assertEquals(expression(), Complex(9.0, 0.0))
} }
} }

@ -1,55 +0,0 @@
# The Core Module (`kmath-core`)
The core features of KMath:
- [algebras](src/commonMain/kotlin/kscience/kmath/operations/Algebra.kt) : Algebraic structures like rings, spaces and fields.
- [nd](src/commonMain/kotlin/kscience/kmath/structures/NDStructure.kt) : Many-dimensional structures and operations on them.
- [linear](src/commonMain/kotlin/kscience/kmath/operations/Algebra.kt) : Basic linear algebra operations (sums, products, etc.), backed by the `Space` API. Advanced linear algebra operations like matrix inversion and LU decomposition.
- [buffers](src/commonMain/kotlin/kscience/kmath/structures/Buffers.kt) : One-dimensional structure
- [expressions](src/commonMain/kotlin/kscience/kmath/expressions) : By writing a single mathematical expression once, users will be able to apply different types of
objects to the expression by providing a context. Expressions can be used for a wide variety of purposes from high
performance calculations to code generation.
- [domains](src/commonMain/kotlin/kscience/kmath/domains) : Domains
- [autodif](src/commonMain/kotlin/kscience/kmath/expressions/SimpleAutoDiff.kt) : Automatic differentiation
> #### Artifact:
>
> This module artifact: `space.kscience:kmath-core:0.2.0`.
>
> Bintray release version: [ ![Download](https://api.bintray.com/packages/mipt-npm/kscience/kmath-core/images/download.svg) ](https://bintray.com/mipt-npm/kscience/kmath-core/_latestVersion)
>
> Bintray development version: [ ![Download](https://api.bintray.com/packages/mipt-npm/dev/kmath-core/images/download.svg) ](https://bintray.com/mipt-npm/dev/kmath-core/_latestVersion)
>
> **Gradle:**
>
> ```gradle
> repositories {
> maven { url 'https://repo.kotlin.link' }
> maven { url 'https://dl.bintray.com/hotkeytlt/maven' }
> maven { url "https://dl.bintray.com/kotlin/kotlin-eap" } // include for builds based on kotlin-eap
>// Uncomment if repo.kotlin.link is unavailable
>// maven { url 'https://dl.bintray.com/mipt-npm/kscience' }
>// maven { url 'https://dl.bintray.com/mipt-npm/dev' }
> }
>
> dependencies {
> implementation 'space.kscience:kmath-core:0.2.0'
> }
> ```
> **Gradle Kotlin DSL:**
>
> ```kotlin
> repositories {
> maven("https://repo.kotlin.link")
> maven("https://dl.bintray.com/kotlin/kotlin-eap") // include for builds based on kotlin-eap
> maven("https://dl.bintray.com/hotkeytlt/maven") // required for a
>// Uncomment if repo.kotlin.link is unavailable
>// maven("https://dl.bintray.com/mipt-npm/kscience")
>// maven("https://dl.bintray.com/mipt-npm/dev")
> }
>
> dependencies {
> implementation("space.kscience:kmath-core:0.2.0")
> }
> ```

File diff suppressed because it is too large Load Diff

@ -16,8 +16,8 @@
package space.kscience.kmath.domains package space.kscience.kmath.domains
import space.kscience.kmath.linear.Point import space.kscience.kmath.linear.Point
import space.kscience.kmath.misc.UnstableKMathAPI
import space.kscience.kmath.structures.Buffer import space.kscience.kmath.structures.Buffer
import space.kscience.kmath.structures.RealBuffer
import space.kscience.kmath.structures.indices import space.kscience.kmath.structures.indices
/** /**
@ -26,6 +26,7 @@ import space.kscience.kmath.structures.indices
* *
* @author Alexander Nozik * @author Alexander Nozik
*/ */
@UnstableKMathAPI
public class HyperSquareDomain(private val lower: Buffer<Double>, private val upper: Buffer<Double>) : RealDomain { public class HyperSquareDomain(private val lower: Buffer<Double>, private val upper: Buffer<Double>) : RealDomain {
public override val dimension: Int get() = lower.size public override val dimension: Int get() = lower.size
@ -33,26 +34,10 @@ public class HyperSquareDomain(private val lower: Buffer<Double>, private val up
point[i] in lower[i]..upper[i] point[i] in lower[i]..upper[i]
} }
public override fun getLowerBound(num: Int, point: Point<Double>): Double = lower[num]
public override fun getLowerBound(num: Int): Double = lower[num] public override fun getLowerBound(num: Int): Double = lower[num]
public override fun getUpperBound(num: Int, point: Point<Double>): Double = upper[num]
public override fun getUpperBound(num: Int): Double = upper[num] public override fun getUpperBound(num: Int): Double = upper[num]
public override fun nearestInDomain(point: Point<Double>): Point<Double> {
val res = DoubleArray(point.size) { i ->
when {
point[i] < lower[i] -> lower[i]
point[i] > upper[i] -> upper[i]
else -> point[i]
}
}
return RealBuffer(*res)
}
public override fun volume(): Double { public override fun volume(): Double {
var res = 1.0 var res = 1.0

@ -15,45 +15,27 @@
*/ */
package space.kscience.kmath.domains package space.kscience.kmath.domains
import space.kscience.kmath.linear.Point import space.kscience.kmath.misc.UnstableKMathAPI
/** /**
* n-dimensional volume * n-dimensional volume
* *
* @author Alexander Nozik * @author Alexander Nozik
*/ */
@UnstableKMathAPI
public interface RealDomain : Domain<Double> { public interface RealDomain : Domain<Double> {
public fun nearestInDomain(point: Point<Double>): Point<Double>
/**
* The lower edge for the domain going down from point
* @param num
* @param point
* @return
*/
public fun getLowerBound(num: Int, point: Point<Double>): Double?
/**
* The upper edge of the domain going up from point
* @param num
* @param point
* @return
*/
public fun getUpperBound(num: Int, point: Point<Double>): Double?
/** /**
* Global lower edge * Global lower edge
* @param num * @param num axis number
* @return
*/ */
public fun getLowerBound(num: Int): Double? public fun getLowerBound(num: Int): Double
/** /**
* Global upper edge * Global upper edge
* @param num * @param num axis number
* @return
*/ */
public fun getUpperBound(num: Int): Double? public fun getUpperBound(num: Int): Double
/** /**
* Hyper volume * Hyper volume

@ -16,19 +16,15 @@
package space.kscience.kmath.domains package space.kscience.kmath.domains
import space.kscience.kmath.linear.Point import space.kscience.kmath.linear.Point
import space.kscience.kmath.misc.UnstableKMathAPI
@UnstableKMathAPI
public class UnconstrainedDomain(public override val dimension: Int) : RealDomain { public class UnconstrainedDomain(public override val dimension: Int) : RealDomain {
public override operator fun contains(point: Point<Double>): Boolean = true public override operator fun contains(point: Point<Double>): Boolean = true
public override fun getLowerBound(num: Int, point: Point<Double>): Double? = Double.NEGATIVE_INFINITY public override fun getLowerBound(num: Int): Double = Double.NEGATIVE_INFINITY
public override fun getLowerBound(num: Int): Double? = Double.NEGATIVE_INFINITY public override fun getUpperBound(num: Int): Double = Double.POSITIVE_INFINITY
public override fun getUpperBound(num: Int, point: Point<Double>): Double? = Double.POSITIVE_INFINITY
public override fun getUpperBound(num: Int): Double? = Double.POSITIVE_INFINITY
public override fun nearestInDomain(point: Point<Double>): Point<Double> = point
public override fun volume(): Double = Double.POSITIVE_INFINITY public override fun volume(): Double = Double.POSITIVE_INFINITY
} }

@ -1,11 +1,11 @@
package space.kscience.kmath.domains package space.kscience.kmath.domains
import space.kscience.kmath.linear.Point import space.kscience.kmath.linear.Point
import space.kscience.kmath.structures.asBuffer import space.kscience.kmath.misc.UnstableKMathAPI
@UnstableKMathAPI
public inline class UnivariateDomain(public val range: ClosedFloatingPointRange<Double>) : RealDomain { public inline class UnivariateDomain(public val range: ClosedFloatingPointRange<Double>) : RealDomain {
public override val dimension: Int public override val dimension: Int get() = 1
get() = 1
public operator fun contains(d: Double): Boolean = range.contains(d) public operator fun contains(d: Double): Boolean = range.contains(d)
@ -14,33 +14,12 @@ public inline class UnivariateDomain(public val range: ClosedFloatingPointRange<
return contains(point[0]) return contains(point[0])
} }
public override fun nearestInDomain(point: Point<Double>): Point<Double> { public override fun getLowerBound(num: Int): Double {
require(point.size == 1)
val value = point[0]
return when {
value in range -> point
value >= range.endInclusive -> doubleArrayOf(range.endInclusive).asBuffer()
else -> doubleArrayOf(range.start).asBuffer()
}
}
public override fun getLowerBound(num: Int, point: Point<Double>): Double? {
require(num == 0) require(num == 0)
return range.start return range.start
} }
public override fun getUpperBound(num: Int, point: Point<Double>): Double? { public override fun getUpperBound(num: Int): Double {
require(num == 0)
return range.endInclusive
}
public override fun getLowerBound(num: Int): Double? {
require(num == 0)
return range.start
}
public override fun getUpperBound(num: Int): Double? {
require(num == 0) require(num == 0)
return range.endInclusive return range.endInclusive
} }

@ -41,25 +41,28 @@ public abstract class FunctionalExpressionAlgebra<T, A : Algebra<T>>(
} }
/** /**
* A context class for [Expression] construction for [Space] algebras. * A context class for [Expression] construction for [Group] algebras.
*/ */
public open class FunctionalExpressionSpace<T, A : Space<T>>( public open class FunctionalExpressionGroup<T, A : Group<T>>(
algebra: A, algebra: A,
) : FunctionalExpressionAlgebra<T, A>(algebra), Space<Expression<T>> { ) : FunctionalExpressionAlgebra<T, A>(algebra), Group<Expression<T>> {
public override val zero: Expression<T> get() = const(algebra.zero) public override val zero: Expression<T> get() = const(algebra.zero)
override fun Expression<T>.unaryMinus(): Expression<T> =
unaryOperation(GroupOperations.MINUS_OPERATION, this)
/** /**
* Builds an Expression of addition of two another expressions. * Builds an Expression of addition of two another expressions.
*/ */
public override fun add(a: Expression<T>, b: Expression<T>): Expression<T> = public override fun add(a: Expression<T>, b: Expression<T>): Expression<T> =
binaryOperationFunction(SpaceOperations.PLUS_OPERATION)(a, b) binaryOperation(GroupOperations.PLUS_OPERATION, a, b)
/** // /**
* Builds an Expression of multiplication of expression by number. // * Builds an Expression of multiplication of expression by number.
*/ // */
public override fun multiply(a: Expression<T>, k: Number): Expression<T> = Expression { arguments -> // public override fun multiply(a: Expression<T>, k: Number): Expression<T> = Expression { arguments ->
algebra.multiply(a.invoke(arguments), k) // algebra.multiply(a.invoke(arguments), k)
} // }
public operator fun Expression<T>.plus(arg: T): Expression<T> = this + const(arg) public operator fun Expression<T>.plus(arg: T): Expression<T> = this + const(arg)
public operator fun Expression<T>.minus(arg: T): Expression<T> = this - const(arg) public operator fun Expression<T>.minus(arg: T): Expression<T> = this - const(arg)
@ -71,13 +74,13 @@ public open class FunctionalExpressionSpace<T, A : Space<T>>(
public override fun binaryOperationFunction(operation: String): (left: Expression<T>, right: Expression<T>) -> Expression<T> = public override fun binaryOperationFunction(operation: String): (left: Expression<T>, right: Expression<T>) -> Expression<T> =
super<FunctionalExpressionAlgebra>.binaryOperationFunction(operation) super<FunctionalExpressionAlgebra>.binaryOperationFunction(operation)
} }
public open class FunctionalExpressionRing<T, A : Ring<T>>( public open class FunctionalExpressionRing<T, A : Ring<T>>(
algebra: A, algebra: A,
) : FunctionalExpressionSpace<T, A>(algebra), Ring<Expression<T>> { ) : FunctionalExpressionGroup<T, A>(algebra), Ring<Expression<T>> {
public override val one: Expression<T> public override val one: Expression<T> get() = const(algebra.one)
get() = const(algebra.one)
/** /**
* Builds an Expression of multiplication of two expressions. * Builds an Expression of multiplication of two expressions.
@ -89,15 +92,16 @@ public open class FunctionalExpressionRing<T, A : Ring<T>>(
public operator fun T.times(arg: Expression<T>): Expression<T> = arg * this public operator fun T.times(arg: Expression<T>): Expression<T> = arg * this
public override fun unaryOperationFunction(operation: String): (arg: Expression<T>) -> Expression<T> = public override fun unaryOperationFunction(operation: String): (arg: Expression<T>) -> Expression<T> =
super<FunctionalExpressionSpace>.unaryOperationFunction(operation) super<FunctionalExpressionGroup>.unaryOperationFunction(operation)
public override fun binaryOperationFunction(operation: String): (left: Expression<T>, right: Expression<T>) -> Expression<T> = public override fun binaryOperationFunction(operation: String): (left: Expression<T>, right: Expression<T>) -> Expression<T> =
super<FunctionalExpressionSpace>.binaryOperationFunction(operation) super<FunctionalExpressionGroup>.binaryOperationFunction(operation)
} }
public open class FunctionalExpressionField<T, A : Field<T>>( public open class FunctionalExpressionField<T, A: Field<T>>(
algebra: A, algebra: A,
) : FunctionalExpressionRing<T, A>(algebra), Field<Expression<T>> { ) : FunctionalExpressionRing<T, A>(algebra), Field<Expression<T>>,
ScaleOperations<Expression<T>> {
/** /**
* Builds an Expression of division an expression by another one. * Builds an Expression of division an expression by another one.
*/ */
@ -112,6 +116,10 @@ public open class FunctionalExpressionField<T, A : Field<T>>(
public override fun binaryOperationFunction(operation: String): (left: Expression<T>, right: Expression<T>) -> Expression<T> = public override fun binaryOperationFunction(operation: String): (left: Expression<T>, right: Expression<T>) -> Expression<T> =
super<FunctionalExpressionRing>.binaryOperationFunction(operation) super<FunctionalExpressionRing>.binaryOperationFunction(operation)
override fun scale(a: Expression<T>, value: Double): Expression<T> = algebra {
Expression { args -> a(args) * value }
}
} }
public open class FunctionalExpressionExtendedField<T, A : ExtendedField<T>>( public open class FunctionalExpressionExtendedField<T, A : ExtendedField<T>>(
@ -151,8 +159,8 @@ public open class FunctionalExpressionExtendedField<T, A : ExtendedField<T>>(
super<FunctionalExpressionField>.binaryOperationFunction(operation) super<FunctionalExpressionField>.binaryOperationFunction(operation)
} }
public inline fun <T, A : Space<T>> A.expressionInSpace(block: FunctionalExpressionSpace<T, A>.() -> Expression<T>): Expression<T> = public inline fun <T, A : Group<T>> A.expressionInSpace(block: FunctionalExpressionGroup<T, A>.() -> Expression<T>): Expression<T> =
FunctionalExpressionSpace(this).block() FunctionalExpressionGroup(this).block()
public inline fun <T, A : Ring<T>> A.expressionInRing(block: FunctionalExpressionRing<T, A>.() -> Expression<T>): Expression<T> = public inline fun <T, A : Ring<T>> A.expressionInRing(block: FunctionalExpressionRing<T, A>.() -> Expression<T>): Expression<T> =
FunctionalExpressionRing(this).block() FunctionalExpressionRing(this).block()
@ -160,5 +168,6 @@ public inline fun <T, A : Ring<T>> A.expressionInRing(block: FunctionalExpressio
public inline fun <T, A : Field<T>> A.expressionInField(block: FunctionalExpressionField<T, A>.() -> Expression<T>): Expression<T> = public inline fun <T, A : Field<T>> A.expressionInField(block: FunctionalExpressionField<T, A>.() -> Expression<T>): Expression<T> =
FunctionalExpressionField(this).block() FunctionalExpressionField(this).block()
public inline fun <T, A : ExtendedField<T>> A.expressionInExtendedField(block: FunctionalExpressionExtendedField<T, A>.() -> Expression<T>): Expression<T> = public inline fun <T, A : ExtendedField<T>> A.expressionInExtendedField(
FunctionalExpressionExtendedField(this).block() block: FunctionalExpressionExtendedField<T, A>.() -> Expression<T>,
): Expression<T> = FunctionalExpressionExtendedField(this).block()

@ -47,36 +47,6 @@ public fun <T : Any> DerivationResult<T>.grad(vararg variables: Symbol): Point<T
return variables.map(::derivative).asBuffer() return variables.map(::derivative).asBuffer()
} }
/**
* Runs differentiation and establishes [SimpleAutoDiffField] context inside the block of code.
*
* The partial derivatives are placed in argument `d` variable
*
* Example:
* ```
* val x by symbol // define variable(s) and their values
* val y = RealField.withAutoDiff() { sqr(x) + 5 * x + 3 } // write formulate in deriv context
* assertEquals(17.0, y.x) // the value of result (y)
* assertEquals(9.0, x.d) // dy/dx
* ```
*
* @param body the action in [SimpleAutoDiffField] context returning [AutoDiffVariable] to differentiate with respect to.
* @return the result of differentiation.
*/
public fun <T : Any, F : Field<T>> F.simpleAutoDiff(
bindings: Map<Symbol, T>,
body: SimpleAutoDiffField<T, F>.() -> AutoDiffValue<T>,
): DerivationResult<T> {
contract { callsInPlace(body, InvocationKind.EXACTLY_ONCE) }
return SimpleAutoDiffField(this, bindings).differentiate(body)
}
public fun <T : Any, F : Field<T>> F.simpleAutoDiff(
vararg bindings: Pair<Symbol, T>,
body: SimpleAutoDiffField<T, F>.() -> AutoDiffValue<T>,
): DerivationResult<T> = simpleAutoDiff(bindings.toMap(), body)
/** /**
* Represents field in context of which functions can be derived. * Represents field in context of which functions can be derived.
*/ */
@ -84,12 +54,9 @@ public fun <T : Any, F : Field<T>> F.simpleAutoDiff(
public open class SimpleAutoDiffField<T : Any, F : Field<T>>( public open class SimpleAutoDiffField<T : Any, F : Field<T>>(
public val context: F, public val context: F,
bindings: Map<Symbol, T>, bindings: Map<Symbol, T>,
) : Field<AutoDiffValue<T>>, ExpressionAlgebra<T, AutoDiffValue<T>>, RingWithNumbers<AutoDiffValue<T>> { ) : Field<AutoDiffValue<T>>, ExpressionAlgebra<T, AutoDiffValue<T>>, NumbersAddOperations<AutoDiffValue<T>> {
public override val zero: AutoDiffValue<T> public override val zero: AutoDiffValue<T> get() = const(context.zero)
get() = const(context.zero) public override val one: AutoDiffValue<T> get() = const(context.one)
public override val one: AutoDiffValue<T>
get() = const(context.one)
// this stack contains pairs of blocks and values to apply them to // this stack contains pairs of blocks and values to apply them to
private var stack: Array<Any?> = arrayOfNulls<Any?>(8) private var stack: Array<Any?> = arrayOfNulls<Any?>(8)
@ -137,6 +104,8 @@ public open class SimpleAutoDiffField<T : Any, F : Field<T>>(
override fun const(value: T): AutoDiffValue<T> = AutoDiffValue(value) override fun const(value: T): AutoDiffValue<T> = AutoDiffValue(value)
override fun number(value: Number): AutoDiffValue<T> = const { one * value }
/** /**
* A variable accessing inner state of derivatives. * A variable accessing inner state of derivatives.
* Use this value in inner builders to avoid creating additional derivative bindings. * Use this value in inner builders to avoid creating additional derivative bindings.
@ -175,21 +144,24 @@ public open class SimpleAutoDiffField<T : Any, F : Field<T>>(
return DerivationResult(result.value, bindings.mapValues { it.value.d }, context) return DerivationResult(result.value, bindings.mapValues { it.value.d }, context)
} }
// Overloads for Double constants // // Overloads for Double constants
//
// public override operator fun Number.plus(b: AutoDiffValue<T>): AutoDiffValue<T> =
// derive(const { this@plus.toDouble() * one + b.value }) { z ->
// b.d += z.d
// }
//
// public override operator fun AutoDiffValue<T>.plus(b: Number): AutoDiffValue<T> = b.plus(this)
//
// public override operator fun Number.minus(b: AutoDiffValue<T>): AutoDiffValue<T> =
// derive(const { this@minus.toDouble() * one - b.value }) { z -> b.d -= z.d }
//
// public override operator fun AutoDiffValue<T>.minus(b: Number): AutoDiffValue<T> =
// derive(const { this@minus.value - one * b.toDouble() }) { z -> d += z.d }
public override operator fun Number.plus(b: AutoDiffValue<T>): AutoDiffValue<T> =
derive(const { this@plus.toDouble() * one + b.value }) { z ->
b.d += z.d
}
public override operator fun AutoDiffValue<T>.plus(b: Number): AutoDiffValue<T> = b.plus(this)
public override operator fun Number.minus(b: AutoDiffValue<T>): AutoDiffValue<T> =
derive(const { this@minus.toDouble() * one - b.value }) { z -> b.d -= z.d }
public override operator fun AutoDiffValue<T>.minus(b: Number): AutoDiffValue<T> =
derive(const { this@minus.value - one * b.toDouble() }) { z -> this@minus.d += z.d }
override fun AutoDiffValue<T>.unaryMinus(): AutoDiffValue<T> =
derive(const { -value }) { z -> d -= z.d }
// Basic math (+, -, *, /) // Basic math (+, -, *, /)
@ -211,12 +183,44 @@ public open class SimpleAutoDiffField<T : Any, F : Field<T>>(
b.d -= z.d * a.value / (b.value * b.value) b.d -= z.d * a.value / (b.value * b.value)
} }
public override fun multiply(a: AutoDiffValue<T>, k: Number): AutoDiffValue<T> = public override fun scale(a: AutoDiffValue<T>, value: Double): AutoDiffValue<T> =
derive(const { k.toDouble() * a.value }) { z -> derive(const { value * a.value }) { z ->
a.d += z.d * k.toDouble() a.d += z.d * value
} }
} }
/**
* Runs differentiation and establishes [SimpleAutoDiffField] context inside the block of code.
*
* The partial derivatives are placed in argument `d` variable
*
* Example:
* ```
* val x by symbol // define variable(s) and their values
* val y = RealField.withAutoDiff() { sqr(x) + 5 * x + 3 } // write formulate in deriv context
* assertEquals(17.0, y.x) // the value of result (y)
* assertEquals(9.0, x.d) // dy/dx
* ```
*
* @param body the action in [SimpleAutoDiffField] context returning [AutoDiffVariable] to differentiate with respect to.
* @return the result of differentiation.
*/
public fun <T : Any, F : Field<T>> F.simpleAutoDiff(
bindings: Map<Symbol, T>,
body: SimpleAutoDiffField<T, F>.() -> AutoDiffValue<T>,
): DerivationResult<T> {
contract { callsInPlace(body, InvocationKind.EXACTLY_ONCE) }
return SimpleAutoDiffField(this, bindings).differentiate(body)
}
public fun <T : Any, F : Field<T>> F.simpleAutoDiff(
vararg bindings: Pair<Symbol, T>,
body: SimpleAutoDiffField<T, F>.() -> AutoDiffValue<T>,
): DerivationResult<T> = simpleAutoDiff(bindings.toMap(), body)
/** /**
* A constructs that creates a derivative structure with required order on-demand * A constructs that creates a derivative structure with required order on-demand
*/ */
@ -247,19 +251,20 @@ public fun <T : Any, F : Field<T>> simpleAutoDiff(field: F): AutoDiffProcessor<T
// Extensions for differentiation of various basic mathematical functions // Extensions for differentiation of various basic mathematical functions
// x ^ 2 // x ^ 2
public fun <T : Any, F : Field<T>> SimpleAutoDiffField<T, F>.sqr(x: AutoDiffValue<T>): AutoDiffValue<T> = public fun <T : Any, F : ExtendedField<T>> SimpleAutoDiffField<T, F>.sqr(x: AutoDiffValue<T>): AutoDiffValue<T> =
derive(const { x.value * x.value }) { z -> x.d += z.d * 2 * x.value } derive(const { x.value * x.value }) { z -> x.d += z.d * 2.0 * x.value }
// x ^ 1/2 // x ^ 1/2
public fun <T : Any, F : ExtendedField<T>> SimpleAutoDiffField<T, F>.sqrt(x: AutoDiffValue<T>): AutoDiffValue<T> = public fun <T : Any, F : ExtendedField<T>> SimpleAutoDiffField<T, F>.sqrt(x: AutoDiffValue<T>): AutoDiffValue<T> =
derive(const { sqrt(x.value) }) { z -> x.d += z.d * 0.5 / z.value } derive(const { sqrt(x.value) }) { z -> x.d += z.d / 2.0 / z.value }
// x ^ y (const) // x ^ y (const)
public fun <T : Any, F : ExtendedField<T>> SimpleAutoDiffField<T, F>.pow( public fun <T : Any, F : ExtendedField<T>> SimpleAutoDiffField<T, F>.pow(
x: AutoDiffValue<T>, x: AutoDiffValue<T>,
y: Double, y: Double,
): AutoDiffValue<T> = ): AutoDiffValue<T> = derive(const { power(x.value, y) }) { z ->
derive(const { power(x.value, y) }) { z -> x.d += z.d * y * power(x.value, y - 1) } x.d += z.d * y * power(x.value, y - 1)
}
public fun <T : Any, F : ExtendedField<T>> SimpleAutoDiffField<T, F>.pow( public fun <T : Any, F : ExtendedField<T>> SimpleAutoDiffField<T, F>.pow(
x: AutoDiffValue<T>, x: AutoDiffValue<T>,
@ -328,7 +333,13 @@ public fun <T : Any, F : ExtendedField<T>> SimpleAutoDiffField<T, F>.atanh(x: Au
public class SimpleAutoDiffExtendedField<T : Any, F : ExtendedField<T>>( public class SimpleAutoDiffExtendedField<T : Any, F : ExtendedField<T>>(
context: F, context: F,
bindings: Map<Symbol, T>, bindings: Map<Symbol, T>,
) : ExtendedField<AutoDiffValue<T>>, SimpleAutoDiffField<T, F>(context, bindings) { ) : ExtendedField<AutoDiffValue<T>>, ScaleOperations<AutoDiffValue<T>>,
SimpleAutoDiffField<T, F>(context, bindings) {
override fun number(value: Number): AutoDiffValue<T> = const { number(value) }
override fun scale(a: AutoDiffValue<T>, value: Double): AutoDiffValue<T> = a * number(value)
// x ^ 2 // x ^ 2
public fun sqr(x: AutoDiffValue<T>): AutoDiffValue<T> = public fun sqr(x: AutoDiffValue<T>): AutoDiffValue<T> =
(this as SimpleAutoDiffField<T, F>).sqr(x) (this as SimpleAutoDiffField<T, F>).sqr(x)

@ -2,18 +2,18 @@ package space.kscience.kmath.expressions
import space.kscience.kmath.operations.ExtendedField import space.kscience.kmath.operations.ExtendedField
import space.kscience.kmath.operations.Field import space.kscience.kmath.operations.Field
import space.kscience.kmath.operations.Group
import space.kscience.kmath.operations.Ring import space.kscience.kmath.operations.Ring
import space.kscience.kmath.operations.Space
import kotlin.contracts.InvocationKind import kotlin.contracts.InvocationKind
import kotlin.contracts.contract import kotlin.contracts.contract
/** /**
* Creates a functional expression with this [Space]. * Creates a functional expression with this [Group].
*/ */
public inline fun <T> Space<T>.spaceExpression(block: FunctionalExpressionSpace<T, Space<T>>.() -> Expression<T>): Expression<T> { public inline fun <T> Group<T>.spaceExpression(block: FunctionalExpressionGroup<T, Group<T>>.() -> Expression<T>): Expression<T> {
contract { callsInPlace(block, InvocationKind.EXACTLY_ONCE) } contract { callsInPlace(block, InvocationKind.EXACTLY_ONCE) }
return FunctionalExpressionSpace(this).block() return FunctionalExpressionGroup(this).block()
} }
/** /**

@ -0,0 +1,83 @@
package space.kscience.kmath.linear
import space.kscience.kmath.nd.*
import space.kscience.kmath.operations.Ring
import space.kscience.kmath.operations.invoke
import space.kscience.kmath.structures.Buffer
import space.kscience.kmath.structures.BufferFactory
import space.kscience.kmath.structures.VirtualBuffer
import space.kscience.kmath.structures.indices
public class BufferLinearSpace<T : Any, A : Ring<T>>(
override val elementAlgebra: A,
private val bufferFactory: BufferFactory<T>,
) : LinearSpace<T, A> {
private fun ndRing(
rows: Int,
cols: Int,
): BufferedNDRing<T, A> = NDAlgebra.ring(elementAlgebra, bufferFactory, rows, cols)
override fun buildMatrix(rows: Int, columns: Int, initializer: A.(i: Int, j: Int) -> T): Matrix<T> =
ndRing(rows, columns).produce { (i, j) -> elementAlgebra.initializer(i, j) }.as2D()
override fun buildVector(size: Int, initializer: A.(Int) -> T): Point<T> =
bufferFactory(size) { elementAlgebra.initializer(it) }
override fun Matrix<T>.unaryMinus(): Matrix<T> = ndRing(rowNum, colNum).run {
unwrap().map { -it }.as2D()
}
override fun Matrix<T>.plus(other: Matrix<T>): Matrix<T> = ndRing(rowNum, colNum).run {
require(shape.contentEquals(other.shape)) { "Shape mismatch on Matrix::plus. Expected $shape but found ${other.shape}" }
unwrap().plus(other.unwrap()).as2D()
}
override fun Matrix<T>.minus(other: Matrix<T>): Matrix<T> = ndRing(rowNum, colNum).run {
require(shape.contentEquals(other.shape)) { "Shape mismatch on Matrix::minus. Expected $shape but found ${other.shape}" }
unwrap().minus(other.unwrap()).as2D()
}
private fun Buffer<T>.linearize() = if (this is VirtualBuffer) {
buildVector(size) { get(it) }
} else {
this
}
override fun Matrix<T>.dot(other: Matrix<T>): Matrix<T> {
require(colNum == other.rowNum) { "Matrix dot operation dimension mismatch: ($rowNum, $colNum) x (${other.rowNum}, ${other.colNum})" }
return elementAlgebra {
val rows = this@dot.rows.map{it.linearize()}
val columns = other.columns.map { it.linearize() }
buildMatrix(rowNum, other.colNum) { i, j ->
val r = rows[i]
val c = columns[j]
var res = zero
for (l in r.indices) {
res += r[l] * c[l]
}
res
}
}
}
override fun Matrix<T>.dot(vector: Point<T>): Point<T> {
require(colNum == vector.size) { "Matrix dot vector operation dimension mismatch: ($rowNum, $colNum) x (${vector.size})" }
return elementAlgebra {
val rows = this@dot.rows.map { it.linearize() }
buildVector(rowNum) { i ->
val r = rows[i]
var res = zero
for (j in r.indices) {
res += r[j] * vector[j]
}
res
}
}
}
override fun Matrix<T>.times(value: T): Matrix<T> = ndRing(rowNum, colNum).run {
unwrap().map { it * value }.as2D()
}
}

@ -1,136 +0,0 @@
package space.kscience.kmath.linear
import space.kscience.kmath.nd.NDStructure
import space.kscience.kmath.nd.Structure2D
import space.kscience.kmath.operations.Ring
import space.kscience.kmath.operations.invoke
import space.kscience.kmath.structures.Buffer
import space.kscience.kmath.structures.BufferFactory
import space.kscience.kmath.structures.asSequence
/**
* Alias for [Structure2D] with more familiar name.
*
* @param T the type of items.
*/
public typealias Matrix<T> = Structure2D<T>
/**
* Basic implementation of Matrix space based on [NDStructure]
*/
public class BufferMatrixContext<T : Any, R : Ring<T>>(
public override val elementContext: R,
private val bufferFactory: BufferFactory<T>,
) : GenericMatrixContext<T, R, BufferMatrix<T>> {
public 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)
}
public override fun point(size: Int, initializer: (Int) -> T): Point<T> = bufferFactory(size, initializer)
private fun Matrix<T>.toBufferMatrix(): BufferMatrix<T> = if (this is BufferMatrix) this else {
produce(rowNum, colNum) { i, j -> get(i, j) }
}
public fun one(rows: Int, columns: Int): Matrix<Double> = VirtualMatrix(rows, columns) { i, j ->
if (i == j) 1.0 else 0.0
} + DiagonalFeature
public override infix fun Matrix<T>.dot(other: Matrix<T>): BufferMatrix<T> {
require(colNum == other.rowNum) { "Matrix dot operation dimension mismatch: ($rowNum, $colNum) x (${other.rowNum}, ${other.colNum})" }
val bufferMatrix = toBufferMatrix()
val otherBufferMatrix = other.toBufferMatrix()
return elementContext {
produce(rowNum, other.colNum) { i, j ->
var res = one
for (l in 0 until colNum) {
res += bufferMatrix[i, l] * otherBufferMatrix[l, j]
}
res
}
}
}
public override infix fun Matrix<T>.dot(vector: Point<T>): Point<T> {
require(colNum == vector.size) { "Matrix dot vector operation dimension mismatch: ($rowNum, $colNum) x (${vector.size})" }
val bufferMatrix = toBufferMatrix()
return elementContext {
bufferFactory(rowNum) { i ->
var res = one
for (j in 0 until colNum) {
res += bufferMatrix[i, j] * vector[j]
}
res
}
}
}
override fun add(a: Matrix<T>, b: Matrix<T>): BufferMatrix<T> {
require(a.rowNum == b.rowNum) { "Row number mismatch in matrix addition. Left side: ${a.rowNum}, right side: ${b.rowNum}" }
require(a.colNum == b.colNum) { "Column number mismatch in matrix addition. Left side: ${a.colNum}, right side: ${b.colNum}" }
val aBufferMatrix = a.toBufferMatrix()
val bBufferMatrix = b.toBufferMatrix()
return elementContext {
produce(a.rowNum, a.colNum) { i, j ->
aBufferMatrix[i, j] + bBufferMatrix[i, j]
}
}
}
override fun multiply(a: Matrix<T>, k: Number): BufferMatrix<T> {
val aBufferMatrix = a.toBufferMatrix()
return elementContext {
produce(a.rowNum, a.colNum) { i, j -> aBufferMatrix[i, j] * k.toDouble() }
}
}
public companion object
}
public class BufferMatrix<T : Any>(
public override val rowNum: Int,
public override val colNum: Int,
public val buffer: Buffer<T>,
) : Matrix<T> {
init {
require(buffer.size == rowNum * colNum) { "Dimension mismatch for matrix structure" }
}
override val shape: IntArray get() = intArrayOf(rowNum, colNum)
public override operator fun get(index: IntArray): T = get(index[0], index[1])
public override operator fun get(i: Int, j: Int): T = buffer[i * colNum + j]
public 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))
}
public override fun equals(other: Any?): Boolean {
if (this === other) return true
return when (other) {
is NDStructure<*> -> NDStructure.contentEquals(this, other)
else -> false
}
}
override fun hashCode(): Int {
var result = rowNum
result = 31 * result + colNum
result = 31 * result + buffer.hashCode()
return result
}
public override fun toString(): String {
return if (rowNum <= 5 && colNum <= 5)
"Matrix(rowsNum = $rowNum, colNum = $colNum)\n" +
rows.asSequence().joinToString(prefix = "(", postfix = ")", separator = "\n ") { buffer ->
buffer.asSequence().joinToString(separator = "\t") { it.toString() }
}
else "Matrix(rowsNum = $rowNum, colNum = $colNum)"
}
}

@ -1,25 +1,22 @@
package space.kscience.kmath.linear package space.kscience.kmath.linear
import space.kscience.kmath.structures.Buffer import space.kscience.kmath.nd.as1D
import space.kscience.kmath.structures.VirtualBuffer
public typealias Point<T> = Buffer<T>
/** /**
* A group of methods to resolve equation A dot X = B, where A and B are matrices or vectors * A group of methods to resolve equation A dot X = B, where A and B are matrices or vectors
*/ */
public interface LinearSolver<T : Any> { public interface LinearSolver<T : Any> {
public fun solve(a: Matrix<T>, b: Matrix<T>): Matrix<T> public fun solve(a: Matrix<T>, b: Matrix<T>): Matrix<T>
public fun solve(a: Matrix<T>, b: Point<T>): Point<T> = solve(a, b.asMatrix()).asPoint() public fun solve(a: Matrix<T>, b: Point<T>): Point<T> = solve(a, b.asMatrix()).asVector()
public fun inverse(a: Matrix<T>): Matrix<T> public fun inverse(a: Matrix<T>): Matrix<T>
} }
/** /**
* Convert matrix to vector if it is possible * Convert matrix to vector if it is possible
*/ */
public fun <T : Any> Matrix<T>.asPoint(): Point<T> = public fun <T : Any> Matrix<T>.asVector(): Point<T> =
if (this.colNum == 1) if (this.colNum == 1)
VirtualBuffer(rowNum) { get(it, 0) } as1D()
else else
error("Can't convert matrix with more than one column to vector") error("Can't convert matrix with more than one column to vector")

@ -0,0 +1,204 @@
package space.kscience.kmath.linear
import space.kscience.kmath.misc.UnstableKMathAPI
import space.kscience.kmath.nd.*
import space.kscience.kmath.operations.*
import space.kscience.kmath.structures.Buffer
import space.kscience.kmath.structures.BufferFactory
import kotlin.reflect.KClass
/**
* Alias for [Structure2D] with more familiar name.
*
* @param T the type of items.
*/
public typealias Matrix<T> = Structure2D<T>
/**
* Alias or using [Buffer] as a point/vector in a many-dimensional space.
*/
public typealias Point<T> = Buffer<T>
/**
* Basic operations on matrices and vectors. Operates on [Matrix].
*
* @param T the type of items in the matrices.
* @param M the type of operated matrices.
*/
public interface LinearSpace<T : Any, out A : Ring<T>> {
public val elementAlgebra: A
/**
* Produces a matrix with this context and given dimensions.
*/
public fun buildMatrix(rows: Int, columns: Int, initializer: A.(i: Int, j: Int) -> T): Matrix<T>
/**
* Produces a point compatible with matrix space (and possibly optimized for it).
*/
public fun buildVector(size: Int, initializer: A.(Int) -> T): Point<T>
public operator fun Matrix<T>.unaryMinus(): Matrix<T> = buildMatrix(rowNum, colNum) { i, j ->
-get(i, j)
}
public operator fun Point<T>.unaryMinus(): Point<T> = buildVector(size) {
-get(it)
}
/**
* Matrix sum
*/
public operator fun Matrix<T>.plus(other: Matrix<T>): Matrix<T> = buildMatrix(rowNum, colNum) { i, j ->
get(i, j) + other[i, j]
}
/**
* Vector sum
*/
public operator fun Point<T>.plus(other: Point<T>): Point<T> = buildVector(size) {
get(it) + other[it]
}
/**
* Matrix subtraction
*/
public operator fun Matrix<T>.minus(other: Matrix<T>): Matrix<T> = buildMatrix(rowNum, colNum) { i, j ->
get(i, j) - other[i, j]
}
/**
* Vector subtraction
*/
public operator fun Point<T>.minus(other: Point<T>): Point<T> = buildVector(size) {
get(it) - other[it]
}
/**
* Computes the dot product of this matrix and another one.
*
* @receiver the multiplicand.
* @param other the multiplier.
* @return the dot product.
*/
public infix fun Matrix<T>.dot(other: Matrix<T>): Matrix<T> {
require(colNum == other.rowNum) { "Matrix dot operation dimension mismatch: ($rowNum, $colNum) x (${other.rowNum}, ${other.colNum})" }
return elementAlgebra {
buildMatrix(rowNum, other.colNum) { i, j ->
var res = zero
for (l in 0 until colNum) {
res += this@dot[i, l] * other[l, j]
}
res
}
}
}
/**
* Computes the dot product of this matrix and a vector.
*
* @receiver the multiplicand.
* @param vector the multiplier.
* @return the dot product.
*/
public infix fun Matrix<T>.dot(vector: Point<T>): Point<T> {
require(colNum == vector.size) { "Matrix dot vector operation dimension mismatch: ($rowNum, $colNum) x (${vector.size})" }
return elementAlgebra {
buildVector(rowNum) { i ->
var res = one
for (j in 0 until colNum) {
res += this@dot[i, j] * vector[j]
}
res
}
}
}
/**
* Multiplies a matrix by its element.
*
* @receiver the multiplicand.
* @param value the multiplier.
* @receiver the product.
*/
public operator fun Matrix<T>.times(value: T): Matrix<T> =
buildMatrix(rowNum, colNum) { i, j -> get(i, j) * value }
/**
* Multiplies an element by a matrix of it.
*
* @receiver the multiplicand.
* @param m the multiplier.
* @receiver the product.
*/
public operator fun T.times(m: Matrix<T>): Matrix<T> = m * this
/**
* Multiplies a vector by its element.
*
* @receiver the multiplicand.
* @param value the multiplier.
* @receiver the product.
*/
public operator fun Point<T>.times(value: T): Point<T> =
buildVector(size) { i -> get(i) * value }
/**
* Multiplies an element by a vector of it.
*
* @receiver the multiplicand.
* @param v the multiplier.
* @receiver the product.
*/
public operator fun T.times(v: Point<T>): Point<T> = v * this
/**
* Gets a feature from the matrix. This function may return some additional features to
* [space.kscience.kmath.nd.NDStructure.getFeature].
*
* @param F the type of feature.
* @param m the matrix.
* @param type the [KClass] instance of [F].
* @return a feature object or `null` if it isn't present.
*/
@UnstableKMathAPI
public fun <F : Any> getFeature(m: Matrix<T>, type: KClass<F>): F? = m.getFeature(type)
public companion object {
/**
* A structured matrix with custom buffer
*/
public fun <T : Any, A : Ring<T>> buffered(
algebra: A,
bufferFactory: BufferFactory<T> = Buffer.Companion::boxing,
): LinearSpace<T, A> = BufferLinearSpace(algebra,bufferFactory)
public val real: LinearSpace<Double, RealField> = buffered(RealField, Buffer.Companion::real)
/**
* Automatic buffered matrix, unboxed if it is possible
*/
public inline fun <reified T : Any, A : Ring<T>> auto(ring: A): LinearSpace<T, A> =
buffered(ring, Buffer.Companion::auto)
}
}
public operator fun <LS : LinearSpace<*, *>, R> LS.invoke(block: LS.() -> R): R = run(block)
/**
* Gets a feature from the matrix. This function may return some additional features to
* [space.kscience.kmath.nd.NDStructure.getFeature].
*
* @param T the type of items in the matrices.
* @param M the type of operated matrices.
* @param F the type of feature.
* @receiver the [LinearSpace] of [T].
* @param m the matrix.
* @return a feature object or `null` if it isn't present.
*/
@UnstableKMathAPI
public inline fun <T : Any, reified F : Any> LinearSpace<T, *>.getFeature(m: Matrix<T>): F? = getFeature(m, F::class)

@ -12,7 +12,7 @@ import space.kscience.kmath.structures.MutableBufferFactory
* Common implementation of [LupDecompositionFeature]. * Common implementation of [LupDecompositionFeature].
*/ */
public class LupDecomposition<T : Any>( public class LupDecomposition<T : Any>(
public val context: MatrixContext<T, Matrix<T>>, public val context: LinearSpace<T, *>,
public val elementContext: Field<T>, public val elementContext: Field<T>,
public val lu: Matrix<T>, public val lu: Matrix<T>,
public val pivot: IntArray, public val pivot: IntArray,
@ -62,15 +62,14 @@ public class LupDecomposition<T : Any>(
} }
@PublishedApi @PublishedApi
internal fun <T : Comparable<T>, F : Field<T>> GenericMatrixContext<T, F, *>.abs(value: T): T = internal fun <T : Comparable<T>> LinearSpace<T, Ring<T>>.abs(value: T): T =
if (value > elementContext.zero) value else elementContext { -value } if (value > elementAlgebra.zero) value else elementAlgebra { -value }
/** /**
* Create a lup decomposition of generic matrix. * Create a lup decomposition of generic matrix.
*/ */
public fun <T : Comparable<T>> MatrixContext<T, Matrix<T>>.lup( public fun <T : Comparable<T>> LinearSpace<T, Field<T>>.lup(
factory: MutableBufferFactory<T>, factory: MutableBufferFactory<T>,
elementContext: Field<T>,
matrix: Matrix<T>, matrix: Matrix<T>,
checkSingular: (T) -> Boolean, checkSingular: (T) -> Boolean,
): LupDecomposition<T> { ): LupDecomposition<T> {
@ -80,7 +79,7 @@ public fun <T : Comparable<T>> MatrixContext<T, Matrix<T>>.lup(
//TODO just waits for KEEP-176 //TODO just waits for KEEP-176
BufferAccessor2D(matrix.rowNum, matrix.colNum, factory).run { BufferAccessor2D(matrix.rowNum, matrix.colNum, factory).run {
elementContext { elementAlgebra {
val lu = create(matrix) val lu = create(matrix)
// Initialize permutation array and parity // Initialize permutation array and parity
@ -142,18 +141,18 @@ public fun <T : Comparable<T>> MatrixContext<T, Matrix<T>>.lup(
for (row in col + 1 until m) lu[row, col] /= luDiag for (row in col + 1 until m) lu[row, col] /= luDiag
} }
return LupDecomposition(this@lup, elementContext, lu.collect(), pivot, even) return LupDecomposition(this@lup, elementAlgebra, lu.collect(), pivot, even)
} }
} }
} }
public inline fun <reified T : Comparable<T>, F : Field<T>> GenericMatrixContext<T, F, Matrix<T>>.lup( public inline fun <reified T : Comparable<T>> LinearSpace<T, Field<T>>.lup(
matrix: Matrix<T>, matrix: Matrix<T>,
noinline checkSingular: (T) -> Boolean, noinline checkSingular: (T) -> Boolean,
): LupDecomposition<T> = lup(MutableBuffer.Companion::auto, elementContext, matrix, checkSingular) ): LupDecomposition<T> = lup(MutableBuffer.Companion::auto, matrix, checkSingular)
public fun MatrixContext<Double, Matrix<Double>>.lup(matrix: Matrix<Double>): LupDecomposition<Double> = public fun LinearSpace<Double, RealField>.lup(matrix: Matrix<Double>): LupDecomposition<Double> =
lup(Buffer.Companion::real, RealField, matrix) { it < 1e-11 } lup(Buffer.Companion::real, matrix) { it < 1e-11 }
public fun <T : Any> LupDecomposition<T>.solveWithLup( public fun <T : Any> LupDecomposition<T>.solveWithLup(
factory: MutableBufferFactory<T>, factory: MutableBufferFactory<T>,
@ -198,7 +197,7 @@ public fun <T : Any> LupDecomposition<T>.solveWithLup(
} }
} }
return context.produce(pivot.size, matrix.colNum) { i, j -> bp[i, j] } return context.buildMatrix(pivot.size, matrix.colNum) { i, j -> bp[i, j] }
} }
} }
} }
@ -210,18 +209,18 @@ public inline fun <reified T : Any> LupDecomposition<T>.solveWithLup(matrix: Mat
* Solves a system of linear equations *ax = b** using LUP decomposition. * Solves a system of linear equations *ax = b** using LUP decomposition.
*/ */
@OptIn(UnstableKMathAPI::class) @OptIn(UnstableKMathAPI::class)
public inline fun <reified T : Comparable<T>, F : Field<T>> GenericMatrixContext<T, F, Matrix<T>>.solveWithLup( public inline fun <reified T : Comparable<T>> LinearSpace<T, Field<T>>.solveWithLup(
a: Matrix<T>, a: Matrix<T>,
b: Matrix<T>, b: Matrix<T>,
noinline bufferFactory: MutableBufferFactory<T> = MutableBuffer.Companion::auto, noinline bufferFactory: MutableBufferFactory<T> = MutableBuffer.Companion::auto,
noinline checkSingular: (T) -> Boolean, noinline checkSingular: (T) -> Boolean,
): Matrix<T> { ): Matrix<T> {
// Use existing decomposition if it is provided by matrix // Use existing decomposition if it is provided by matrix
val decomposition = a.getFeature() ?: lup(bufferFactory, elementContext, a, checkSingular) val decomposition = a.getFeature() ?: lup(bufferFactory, a, checkSingular)
return decomposition.solveWithLup(bufferFactory, b) return decomposition.solveWithLup(bufferFactory, b)
} }
public inline fun <reified T : Comparable<T>, F : Field<T>> GenericMatrixContext<T, F, Matrix<T>>.inverseWithLup( public inline fun <reified T : Comparable<T>> LinearSpace<T, Field<T>>.inverseWithLup(
matrix: Matrix<T>, matrix: Matrix<T>,
noinline bufferFactory: MutableBufferFactory<T> = MutableBuffer.Companion::auto, noinline bufferFactory: MutableBufferFactory<T> = MutableBuffer.Companion::auto,
noinline checkSingular: (T) -> Boolean, noinline checkSingular: (T) -> Boolean,
@ -229,15 +228,15 @@ public inline fun <reified T : Comparable<T>, F : Field<T>> GenericMatrixContext
@OptIn(UnstableKMathAPI::class) @OptIn(UnstableKMathAPI::class)
public fun RealMatrixContext.solveWithLup(a: Matrix<Double>, b: Matrix<Double>): Matrix<Double> { public fun LinearSpace<Double, RealField>.solveWithLup(a: Matrix<Double>, b: Matrix<Double>): Matrix<Double> {
// Use existing decomposition if it is provided by matrix // Use existing decomposition if it is provided by matrix
val bufferFactory: MutableBufferFactory<Double> = MutableBuffer.Companion::real val bufferFactory: MutableBufferFactory<Double> = MutableBuffer.Companion::real
val decomposition: LupDecomposition<Double> = a.getFeature() ?: lup(bufferFactory, RealField, a) { it < 1e-11 } val decomposition: LupDecomposition<Double> = a.getFeature() ?: lup(bufferFactory, a) { it < 1e-11 }
return decomposition.solveWithLup(bufferFactory, b) return decomposition.solveWithLup(bufferFactory, b)
} }
/** /**
* Inverses a square matrix using LUP decomposition. Non square matrix will throw a error. * Inverses a square matrix using LUP decomposition. Non square matrix will throw a error.
*/ */
public fun RealMatrixContext.inverseWithLup(matrix: Matrix<Double>): Matrix<Double> = public fun LinearSpace<Double, RealField>.inverseWithLup(matrix: Matrix<Double>): Matrix<Double> =
solveWithLup(matrix, one(matrix.rowNum, matrix.colNum)) solveWithLup(matrix, one(matrix.rowNum, matrix.colNum))

@ -1,46 +1,43 @@
package space.kscience.kmath.linear package space.kscience.kmath.linear
import space.kscience.kmath.nd.Structure2D import space.kscience.kmath.misc.UnstableKMathAPI
import space.kscience.kmath.structures.Buffer import space.kscience.kmath.operations.Ring
import space.kscience.kmath.structures.BufferFactory
import space.kscience.kmath.structures.asBuffer
public class MatrixBuilder(public val rows: Int, public val columns: Int) { public class MatrixBuilder<T : Any, A : Ring<T>>(
public operator fun <T : Any> invoke(vararg elements: T): Matrix<T> { public val linearSpace: LinearSpace<T, A>,
public val rows: Int,
public val columns: Int,
) {
public operator fun invoke(vararg elements: T): Matrix<T> {
require(rows * columns == elements.size) { "The number of elements ${elements.size} is not equal $rows * $columns" } require(rows * columns == elements.size) { "The number of elements ${elements.size} is not equal $rows * $columns" }
val buffer = elements.asBuffer() return linearSpace.buildMatrix(rows, columns) { i, j -> elements[i * columns + j] }
return BufferMatrix(rows, columns, buffer)
} }
//TODO add specific matrix builder functions like diagonal, etc //TODO add specific matrix builder functions like diagonal, etc
} }
public fun Structure2D.Companion.build(rows: Int, columns: Int): MatrixBuilder = MatrixBuilder(rows, columns) /**
* Create a matrix builder with given number of rows and columns
*/
@UnstableKMathAPI
public fun <T : Any, A : Ring<T>> LinearSpace<T, A>.matrix(rows: Int, columns: Int): MatrixBuilder<T, A> =
MatrixBuilder(this, rows, columns)
public fun <T : Any> Structure2D.Companion.row(vararg values: T): Matrix<T> { @UnstableKMathAPI
val buffer = values.asBuffer() public fun <T : Any> LinearSpace<T, Ring<T>>.vector(vararg elements: T): Point<T> {
return BufferMatrix(1, values.size, buffer) return buildVector(elements.size) { elements[it] }
} }
public inline fun <reified T : Any> Structure2D.Companion.row( public inline fun <T : Any> LinearSpace<T, Ring<T>>.row(
size: Int, size: Int,
factory: BufferFactory<T> = Buffer.Companion::auto, crossinline builder: (Int) -> T,
noinline builder: (Int) -> T, ): Matrix<T> = buildMatrix(1, size) { _, j -> builder(j) }
): Matrix<T> {
val buffer = factory(size, builder)
return BufferMatrix(1, size, buffer)
}
public fun <T : Any> Structure2D.Companion.column(vararg values: T): Matrix<T> { public fun <T : Any> LinearSpace<T, Ring<T>>.row(vararg values: T): Matrix<T> = row(values.size, values::get)
val buffer = values.asBuffer()
return BufferMatrix(values.size, 1, buffer)
}
public inline fun <reified T : Any> Structure2D.Companion.column( public inline fun <T : Any> LinearSpace<T, Ring<T>>.column(
size: Int, size: Int,
factory: BufferFactory<T> = Buffer.Companion::auto, crossinline builder: (Int) -> T,
noinline builder: (Int) -> T, ): Matrix<T> = buildMatrix(size, 1) { i, _ -> builder(i) }
): Matrix<T> {
val buffer = factory(size, builder) public fun <T : Any> LinearSpace<T, Ring<T>>.column(vararg values: T): Matrix<T> = column(values.size, values::get)
return BufferMatrix(size, 1, buffer)
}

@ -1,176 +0,0 @@
package space.kscience.kmath.linear
import space.kscience.kmath.misc.UnstableKMathAPI
import space.kscience.kmath.operations.Ring
import space.kscience.kmath.operations.SpaceOperations
import space.kscience.kmath.operations.invoke
import space.kscience.kmath.operations.sum
import space.kscience.kmath.structures.Buffer
import space.kscience.kmath.structures.BufferFactory
import space.kscience.kmath.structures.asSequence
import kotlin.reflect.KClass
/**
* Basic operations on matrices. Operates on [Matrix].
*
* @param T the type of items in the matrices.
* @param M the type of operated matrices.
*/
public interface MatrixContext<T : Any, out M : Matrix<T>> : SpaceOperations<Matrix<T>> {
/**
* Produces a matrix with this context and given dimensions.
*/
public fun produce(rows: Int, columns: Int, initializer: (i: Int, j: Int) -> T): M
/**
* Produces a point compatible with matrix space (and possibly optimized for it).
*/
public fun point(size: Int, initializer: (Int) -> T): Point<T> = Buffer.boxing(size, initializer)
@Suppress("UNCHECKED_CAST")
public override fun binaryOperationFunction(operation: String): (left: Matrix<T>, right: Matrix<T>) -> M =
when (operation) {
"dot" -> { left, right -> left dot right }
else -> super.binaryOperationFunction(operation) as (Matrix<T>, Matrix<T>) -> M
}
/**
* Computes the dot product of this matrix and another one.
*
* @receiver the multiplicand.
* @param other the multiplier.
* @return the dot product.
*/
public infix fun Matrix<T>.dot(other: Matrix<T>): M
/**
* Computes the dot product of this matrix and a vector.
*
* @receiver the multiplicand.
* @param vector the multiplier.
* @return the dot product.
*/
public infix fun Matrix<T>.dot(vector: Point<T>): Point<T>
/**
* Multiplies a matrix by its element.
*
* @receiver the multiplicand.
* @param value the multiplier.
* @receiver the product.
*/
public operator fun Matrix<T>.times(value: T): M
/**
* Multiplies an element by a matrix of it.
*
* @receiver the multiplicand.
* @param m the multiplier.
* @receiver the product.
*/
public operator fun T.times(m: Matrix<T>): M = m * this
/**
* Gets a feature from the matrix. This function may return some additional features to
* [kscience.kmath.nd.NDStructure.getFeature].
*
* @param F the type of feature.
* @param m the matrix.
* @param type the [KClass] instance of [F].
* @return a feature object or `null` if it isn't present.
*/
@UnstableKMathAPI
public fun <F : Any> getFeature(m: Matrix<T>, type: KClass<F>): F? = m.getFeature(type)
public companion object {
/**
* A structured matrix with custom buffer
*/
public fun <T : Any, R : Ring<T>> buffered(
ring: R,
bufferFactory: BufferFactory<T> = Buffer.Companion::boxing,
): GenericMatrixContext<T, R, BufferMatrix<T>> = BufferMatrixContext(ring, bufferFactory)
/**
* Automatic buffered matrix, unboxed if it is possible
*/
public inline fun <reified T : Any, R : Ring<T>> auto(ring: R): GenericMatrixContext<T, R, BufferMatrix<T>> =
buffered(ring, Buffer.Companion::auto)
}
}
/**
* Gets a feature from the matrix. This function may return some additional features to
* [kscience.kmath.nd.NDStructure.getFeature].
*
* @param T the type of items in the matrices.
* @param M the type of operated matrices.
* @param F the type of feature.
* @receiver the [MatrixContext] of [T].
* @param m the matrix.
* @return a feature object or `null` if it isn't present.
*/
@UnstableKMathAPI
public inline fun <T : Any, reified F : Any> MatrixContext<T, *>.getFeature(m: Matrix<T>): F? =
getFeature(m, F::class)
/**
* Partial implementation of [MatrixContext] for matrices of [Ring].
*
* @param T the type of items in the matrices.
* @param R the type of ring of matrix elements.
* @param M the type of operated matrices.
*/
public interface GenericMatrixContext<T : Any, R : Ring<T>, out M : Matrix<T>> : MatrixContext<T, M> {
/**
* The ring over matrix elements.
*/
public val elementContext: R
public override infix fun Matrix<T>.dot(other: Matrix<T>): M {
//TODO add typed error
require(colNum == other.rowNum) { "Matrix dot operation dimension mismatch: ($rowNum, $colNum) x (${other.rowNum}, ${other.colNum})" }
return produce(rowNum, other.colNum) { i, j ->
val row = rows[i]
val column = other.columns[j]
elementContext { sum(row.asSequence().zip(column.asSequence(), ::multiply)) }
}
}
public override infix fun Matrix<T>.dot(vector: Point<T>): Point<T> {
//TODO add typed error
require(colNum == vector.size) { "Matrix dot vector operation dimension mismatch: ($rowNum, $colNum) x (${vector.size})" }
return point(rowNum) { i ->
val row = rows[i]
elementContext { sum(row.asSequence().zip(vector.asSequence(), ::multiply)) }
}
}
public override operator fun Matrix<T>.unaryMinus(): M =
produce(rowNum, colNum) { i, j -> elementContext { -get(i, j) } }
public override fun add(a: Matrix<T>, b: Matrix<T>): M {
require(a.rowNum == b.rowNum && a.colNum == b.colNum) {
"Matrix operation dimension mismatch. [${a.rowNum},${a.colNum}] + [${b.rowNum},${b.colNum}]"
}
return produce(a.rowNum, a.colNum) { i, j -> elementContext { a[i, j] + b[i, j] } }
}
public override operator fun Matrix<T>.minus(b: Matrix<T>): M {
require(rowNum == b.rowNum && colNum == b.colNum) {
"Matrix operation dimension mismatch. [$rowNum,$colNum] - [${b.rowNum},${b.colNum}]"
}
return produce(rowNum, colNum) { i, j -> elementContext { get(i, j) + b[i, j] } }
}
public override fun multiply(a: Matrix<T>, k: Number): M =
produce(a.rowNum, a.colNum) { i, j -> elementContext { a[i, j] * k } }
public override operator fun Matrix<T>.times(value: T): M =
produce(rowNum, colNum) { i, j -> elementContext { get(i, j) * value } }
}

@ -1,20 +1,16 @@
package space.kscience.kmath.linear package space.kscience.kmath.linear
import space.kscience.kmath.misc.UnstableKMathAPI import space.kscience.kmath.misc.UnstableKMathAPI
import space.kscience.kmath.nd.Structure2D
import space.kscience.kmath.nd.getFeature import space.kscience.kmath.nd.getFeature
import space.kscience.kmath.operations.Ring import space.kscience.kmath.operations.Ring
import space.kscience.kmath.structures.asBuffer
import kotlin.math.sqrt
import kotlin.reflect.KClass import kotlin.reflect.KClass
import kotlin.reflect.safeCast
/** /**
* A [Matrix] that holds [MatrixFeature] objects. * A [Matrix] that holds [MatrixFeature] objects.
* *
* @param T the type of items. * @param T the type of items.
*/ */
public class MatrixWrapper<T : Any> internal constructor( public class MatrixWrapper<T : Any> internal constructor(
public val origin: Matrix<T>, public val origin: Matrix<T>,
public val features: Set<MatrixFeature>, public val features: Set<MatrixFeature>,
) : Matrix<T> by origin { ) : Matrix<T> by origin {
@ -23,7 +19,8 @@ public class MatrixWrapper<T : Any> internal constructor(
* Get the first feature matching given class. Does not guarantee that matrix has only one feature matching the criteria * Get the first feature matching given class. Does not guarantee that matrix has only one feature matching the criteria
*/ */
@UnstableKMathAPI @UnstableKMathAPI
override fun <T : Any> getFeature(type: KClass<T>): T? = type.safeCast(features.find { type.isInstance(it) }) @Suppress("UNCHECKED_CAST")
override fun <T : Any> getFeature(type: KClass<T>): T? = features.singleOrNull { type.isInstance(it) } as? T
?: origin.getFeature(type) ?: origin.getFeature(type)
override fun equals(other: Any?): Boolean = origin == other override fun equals(other: Any?): Boolean = origin == other
@ -60,30 +57,26 @@ public operator fun <T : Any> Matrix<T>.plus(newFeatures: Collection<MatrixFeatu
MatrixWrapper(this, newFeatures.toSet()) MatrixWrapper(this, newFeatures.toSet())
} }
/**
* Build a square matrix from given elements.
*/
public fun <T : Any> Structure2D.Companion.square(vararg elements: T): Matrix<T> {
val size: Int = sqrt(elements.size.toDouble()).toInt()
require(size * size == elements.size) { "The number of elements ${elements.size} is not a full square" }
val buffer = elements.asBuffer()
return BufferMatrix(size, size, buffer)
}
/** /**
* Diagonal matrix of ones. The matrix is virtual no actual matrix is created * Diagonal matrix of ones. The matrix is virtual no actual matrix is created
*/ */
public fun <T : Any, R : Ring<T>> GenericMatrixContext<T, R, *>.one(rows: Int, columns: Int): Matrix<T> = public fun <T : Any> LinearSpace<T, Ring<T>>.one(
VirtualMatrix(rows, columns) { i, j -> rows: Int,
if (i == j) elementContext.one else elementContext.zero columns: Int,
} + UnitFeature ): Matrix<T> = VirtualMatrix(rows, columns) { i, j ->
if (i == j) elementAlgebra.one else elementAlgebra.zero
} + UnitFeature
/** /**
* A virtual matrix of zeroes * A virtual matrix of zeroes
*/ */
public fun <T : Any, R : Ring<T>> GenericMatrixContext<T, R, *>.zero(rows: Int, columns: Int): Matrix<T> = public fun <T : Any> LinearSpace<T, Ring<T>>.zero(
VirtualMatrix(rows, columns) { _, _ -> elementContext.zero } + ZeroFeature rows: Int,
columns: Int,
): Matrix<T> = VirtualMatrix(rows, columns) { _, _ ->
elementAlgebra.zero
} + ZeroFeature
public class TransposedFeature<T : Any>(public val original: Matrix<T>) : MatrixFeature public class TransposedFeature<T : Any>(public val original: Matrix<T>) : MatrixFeature
@ -91,9 +84,7 @@ public class TransposedFeature<T : Any>(public val original: Matrix<T>) : Matrix
* Create a virtual transposed matrix without copying anything. `A.transpose().transpose() === A` * Create a virtual transposed matrix without copying anything. `A.transpose().transpose() === A`
*/ */
@OptIn(UnstableKMathAPI::class) @OptIn(UnstableKMathAPI::class)
public fun <T : Any> Matrix<T>.transpose(): Matrix<T> { public fun <T : Any> Matrix<T>.transpose(): Matrix<T> = getFeature<TransposedFeature<T>>()?.original ?: VirtualMatrix(
return getFeature<TransposedFeature<T>>()?.original ?: VirtualMatrix( colNum,
colNum, rowNum,
rowNum, ) { i, j -> get(j, i) } + TransposedFeature(this)
) { i, j -> get(j, i) } + TransposedFeature(this)
}

@ -1,75 +0,0 @@
package space.kscience.kmath.linear
import space.kscience.kmath.structures.RealBuffer
public object RealMatrixContext : MatrixContext<Double, BufferMatrix<Double>> {
public override fun produce(
rows: Int,
columns: Int,
initializer: (i: Int, j: Int) -> Double,
): BufferMatrix<Double> {
val buffer = RealBuffer(rows * columns) { offset -> initializer(offset / columns, offset % columns) }
return BufferMatrix(rows, columns, buffer)
}
public fun Matrix<Double>.toBufferMatrix(): BufferMatrix<Double> = if (this is BufferMatrix) this else {
produce(rowNum, colNum) { i, j -> get(i, j) }
}
public fun one(rows: Int, columns: Int): Matrix<Double> = VirtualMatrix(rows, columns) { i, j ->
if (i == j) 1.0 else 0.0
} + DiagonalFeature
public override infix fun Matrix<Double>.dot(other: Matrix<Double>): BufferMatrix<Double> {
require(colNum == other.rowNum) { "Matrix dot operation dimension mismatch: ($rowNum, $colNum) x (${other.rowNum}, ${other.colNum})" }
val bufferMatrix = toBufferMatrix()
val otherBufferMatrix = other.toBufferMatrix()
return produce(rowNum, other.colNum) { i, j ->
var res = 0.0
for (l in 0 until colNum) {
res += bufferMatrix[i, l] * otherBufferMatrix[l, j]
}
res
}
}
public override infix fun Matrix<Double>.dot(vector: Point<Double>): Point<Double> {
require(colNum == vector.size) { "Matrix dot vector operation dimension mismatch: ($rowNum, $colNum) x (${vector.size})" }
val bufferMatrix = toBufferMatrix()
return RealBuffer(rowNum) { i ->
var res = 0.0
for (j in 0 until colNum) {
res += bufferMatrix[i, j] * vector[j]
}
res
}
}
override fun add(a: Matrix<Double>, b: Matrix<Double>): BufferMatrix<Double> {
require(a.rowNum == b.rowNum) { "Row number mismatch in matrix addition. Left side: ${a.rowNum}, right side: ${b.rowNum}" }
require(a.colNum == b.colNum) { "Column number mismatch in matrix addition. Left side: ${a.colNum}, right side: ${b.colNum}" }
val aBufferMatrix = a.toBufferMatrix()
val bBufferMatrix = b.toBufferMatrix()
return produce(a.rowNum, a.colNum) { i, j ->
aBufferMatrix[i, j] + bBufferMatrix[i, j]
}
}
override fun Matrix<Double>.times(value: Double): BufferMatrix<Double> {
val bufferMatrix = toBufferMatrix()
return produce(rowNum, colNum) { i, j -> bufferMatrix[i, j] * value }
}
override fun multiply(a: Matrix<Double>, k: Number): BufferMatrix<Double> {
val aBufferMatrix = a.toBufferMatrix()
return produce(a.rowNum, a.colNum) { i, j -> aBufferMatrix[i, j] * k.toDouble() }
}
}
/**
* Partially optimized real-valued matrix
*/
public val MatrixContext.Companion.real: RealMatrixContext get() = RealMatrixContext

@ -1,70 +0,0 @@
package space.kscience.kmath.linear
import space.kscience.kmath.operations.RealField
import space.kscience.kmath.operations.Space
import space.kscience.kmath.operations.invoke
import space.kscience.kmath.structures.Buffer
import space.kscience.kmath.structures.BufferFactory
/**
* A linear space for vectors.
* Could be used on any point-like structure
*/
public interface VectorSpace<T : Any, S : Space<T>> : Space<Point<T>> {
public val size: Int
public val space: S
override val zero: Point<T> get() = produce { space.zero }
public fun produce(initializer: S.(Int) -> T): Point<T>
/**
* Produce a space-element of this vector space for expressions
*/
//fun produceElement(initializer: (Int) -> T): Vector<T, S>
override fun add(a: Point<T>, b: Point<T>): Point<T> = produce { space { a[it] + b[it] } }
override fun multiply(a: Point<T>, k: Number): Point<T> = produce { space { a[it] * k } }
//TODO add basis
public companion object {
private val realSpaceCache: MutableMap<Int, BufferVectorSpace<Double, RealField>> = hashMapOf()
/**
* Non-boxing double vector space
*/
public fun real(size: Int): BufferVectorSpace<Double, RealField> = realSpaceCache.getOrPut(size) {
BufferVectorSpace(
size,
RealField,
Buffer.Companion::auto
)
}
/**
* A structured vector space with custom buffer
*/
public fun <T : Any, S : Space<T>> buffered(
size: Int,
space: S,
bufferFactory: BufferFactory<T> = Buffer.Companion::boxing,
): BufferVectorSpace<T, S> = BufferVectorSpace(size, space, bufferFactory)
/**
* Automatic buffered vector, unboxed if it is possible
*/
public inline fun <reified T : Any, S : Space<T>> auto(size: Int, space: S): VectorSpace<T, S> =
buffered(size, space, Buffer.Companion::auto)
}
}
public class BufferVectorSpace<T : Any, S : Space<T>>(
override val size: Int,
override val space: S,
public val bufferFactory: BufferFactory<T>,
) : VectorSpace<T, S> {
override fun produce(initializer: S.(Int) -> T): Buffer<T> = bufferFactory(size) { space.initializer(it) }
//override fun produceElement(initializer: (Int) -> T): Vector<T, S> = BufferVector(this, produce(initializer))
}

@ -1,5 +1,7 @@
package space.kscience.kmath.linear package space.kscience.kmath.linear
import space.kscience.kmath.nd.NDStructure
/** /**
* The matrix where each element is evaluated each time when is being accessed. * The matrix where each element is evaluated each time when is being accessed.
* *
@ -17,12 +19,8 @@ public class VirtualMatrix<T : Any>(
override fun equals(other: Any?): Boolean { override fun equals(other: Any?): Boolean {
if (this === other) return true if (this === other) return true
if (other !is Matrix<*>) return false if (other !is NDStructure<*>) return false
return NDStructure.contentEquals(this, other)
if (rowNum != other.rowNum) return false
if (colNum != other.colNum) return false
return elements().all { (index, value) -> value == other[index] }
} }
override fun hashCode(): Int { override fun hashCode(): Int {
@ -31,6 +29,4 @@ public class VirtualMatrix<T : Any>(
result = 31 * result + generator.hashCode() result = 31 * result + generator.hashCode()
return result return result
} }
} }

@ -1,6 +1,6 @@
package space.kscience.kmath.misc package space.kscience.kmath.misc
import space.kscience.kmath.operations.Space import space.kscience.kmath.operations.Group
import space.kscience.kmath.operations.invoke import space.kscience.kmath.operations.invoke
import kotlin.jvm.JvmName import kotlin.jvm.JvmName
@ -37,7 +37,7 @@ public fun <T, R> List<T>.cumulative(initial: R, operation: (R, T) -> R): List<R
/** /**
* Cumulative sum with custom space * Cumulative sum with custom space
*/ */
public fun <T> Iterable<T>.cumulativeSum(space: Space<T>): Iterable<T> = public fun <T> Iterable<T>.cumulativeSum(space: Group<T>): Iterable<T> =
space { cumulative(zero) { element: T, sum: T -> sum + element } } space { cumulative(zero) { element: T, sum: T -> sum + element } }
@JvmName("cumulativeSumOfDouble") @JvmName("cumulativeSumOfDouble")
@ -49,7 +49,7 @@ public fun Iterable<Int>.cumulativeSum(): Iterable<Int> = cumulative(0) { elemen
@JvmName("cumulativeSumOfLong") @JvmName("cumulativeSumOfLong")
public fun Iterable<Long>.cumulativeSum(): Iterable<Long> = cumulative(0L) { element, sum -> sum + element } public fun Iterable<Long>.cumulativeSum(): Iterable<Long> = cumulative(0L) { element, sum -> sum + element }
public fun <T> Sequence<T>.cumulativeSum(space: Space<T>): Sequence<T> = public fun <T> Sequence<T>.cumulativeSum(space: Group<T>): Sequence<T> =
space { cumulative(zero) { element: T, sum: T -> sum + element } } space { cumulative(zero) { element: T, sum: T -> sum + element } }
@JvmName("cumulativeSumOfDouble") @JvmName("cumulativeSumOfDouble")
@ -61,7 +61,7 @@ public fun Sequence<Int>.cumulativeSum(): Sequence<Int> = cumulative(0) { elemen
@JvmName("cumulativeSumOfLong") @JvmName("cumulativeSumOfLong")
public fun Sequence<Long>.cumulativeSum(): Sequence<Long> = cumulative(0L) { element, sum -> sum + element } public fun Sequence<Long>.cumulativeSum(): Sequence<Long> = cumulative(0L) { element, sum -> sum + element }
public fun <T> List<T>.cumulativeSum(space: Space<T>): List<T> = public fun <T> List<T>.cumulativeSum(space: Group<T>): List<T> =
space { cumulative(zero) { element: T, sum: T -> sum + element } } space { cumulative(zero) { element: T, sum: T -> sum + element } }
@JvmName("cumulativeSumOfDouble") @JvmName("cumulativeSumOfDouble")

@ -1,19 +1,16 @@
package space.kscience.kmath.nd package space.kscience.kmath.nd
import space.kscience.kmath.operations.Field import space.kscience.kmath.operations.*
import space.kscience.kmath.operations.RealField
import space.kscience.kmath.operations.Ring
import space.kscience.kmath.operations.Space
import space.kscience.kmath.structures.Buffer import space.kscience.kmath.structures.Buffer
import space.kscience.kmath.structures.BufferFactory import space.kscience.kmath.structures.BufferFactory
import kotlin.contracts.InvocationKind import kotlin.contracts.InvocationKind
import kotlin.contracts.contract import kotlin.contracts.contract
public interface BufferNDAlgebra<T, C> : NDAlgebra<T, C> { public interface BufferNDAlgebra<T, A : Algebra<T>> : NDAlgebra<T, A> {
public val strides: Strides public val strides: Strides
public val bufferFactory: BufferFactory<T> public val bufferFactory: BufferFactory<T>
override fun produce(initializer: C.(IntArray) -> T): NDBuffer<T> = NDBuffer( override fun produce(initializer: A.(IntArray) -> T): NDBuffer<T> = NDBuffer(
strides, strides,
bufferFactory(strides.linearSize) { offset -> bufferFactory(strides.linearSize) { offset ->
elementContext.initializer(strides.index(offset)) elementContext.initializer(strides.index(offset))
@ -30,14 +27,14 @@ public interface BufferNDAlgebra<T, C> : NDAlgebra<T, C> {
else -> bufferFactory(strides.linearSize) { offset -> get(strides.index(offset)) } else -> bufferFactory(strides.linearSize) { offset -> get(strides.index(offset)) }
} }
override fun NDStructure<T>.map(transform: C.(T) -> T): NDBuffer<T> { override fun NDStructure<T>.map(transform: A.(T) -> T): NDBuffer<T> {
val buffer = bufferFactory(strides.linearSize) { offset -> val buffer = bufferFactory(strides.linearSize) { offset ->
elementContext.transform(buffer[offset]) elementContext.transform(buffer[offset])
} }
return NDBuffer(strides, buffer) return NDBuffer(strides, buffer)
} }
override fun NDStructure<T>.mapIndexed(transform: C.(index: IntArray, T) -> T): NDBuffer<T> { override fun NDStructure<T>.mapIndexed(transform: A.(index: IntArray, T) -> T): NDBuffer<T> {
val buffer = bufferFactory(strides.linearSize) { offset -> val buffer = bufferFactory(strides.linearSize) { offset ->
elementContext.transform( elementContext.transform(
strides.index(offset), strides.index(offset),
@ -47,7 +44,7 @@ public interface BufferNDAlgebra<T, C> : NDAlgebra<T, C> {
return NDBuffer(strides, buffer) return NDBuffer(strides, buffer)
} }
override fun combine(a: NDStructure<T>, b: NDStructure<T>, transform: C.(T, T) -> T): NDBuffer<T> { override fun combine(a: NDStructure<T>, b: NDStructure<T>, transform: A.(T, T) -> T): NDBuffer<T> {
val buffer = bufferFactory(strides.linearSize) { offset -> val buffer = bufferFactory(strides.linearSize) { offset ->
elementContext.transform(a.buffer[offset], b.buffer[offset]) elementContext.transform(a.buffer[offset], b.buffer[offset])
} }
@ -55,20 +52,21 @@ public interface BufferNDAlgebra<T, C> : NDAlgebra<T, C> {
} }
} }
public open class BufferedNDSpace<T, R : Space<T>>( public open class BufferedNDGroup<T, A : Group<T>>(
final override val shape: IntArray, final override val shape: IntArray,
final override val elementContext: R, final override val elementContext: A,
final override val bufferFactory: BufferFactory<T>, final override val bufferFactory: BufferFactory<T>,
) : NDSpace<T, R>, BufferNDAlgebra<T, R> { ) : NDGroup<T, A>, BufferNDAlgebra<T, A> {
override val strides: Strides = DefaultStrides(shape) override val strides: Strides = DefaultStrides(shape)
override val zero: NDBuffer<T> by lazy { produce { zero } } override val zero: NDBuffer<T> by lazy { produce { zero } }
override fun NDStructure<T>.unaryMinus(): NDStructure<T> = produce { -get(it) }
} }
public open class BufferedNDRing<T, R : Ring<T>>( public open class BufferedNDRing<T, R : Ring<T>>(
shape: IntArray, shape: IntArray,
elementContext: R, elementContext: R,
bufferFactory: BufferFactory<T>, bufferFactory: BufferFactory<T>,
) : BufferedNDSpace<T, R>(shape, elementContext, bufferFactory), NDRing<T, R> { ) : BufferedNDGroup<T, R>(shape, elementContext, bufferFactory), NDRing<T, R> {
override val one: NDBuffer<T> by lazy { produce { one } } override val one: NDBuffer<T> by lazy { produce { one } }
} }
@ -76,22 +74,25 @@ public open class BufferedNDField<T, R : Field<T>>(
shape: IntArray, shape: IntArray,
elementContext: R, elementContext: R,
bufferFactory: BufferFactory<T>, bufferFactory: BufferFactory<T>,
) : BufferedNDRing<T, R>(shape, elementContext, bufferFactory), NDField<T, R> ) : BufferedNDRing<T, R>(shape, elementContext, bufferFactory), NDField<T, R> {
// space factories override fun scale(a: NDStructure<T>, value: Double): NDStructure<T> = a.map { it * value }
public fun <T, A : Space<T>> NDAlgebra.Companion.space( }
// group factories
public fun <T, A : Group<T>> NDAlgebra.Companion.group(
space: A, space: A,
bufferFactory: BufferFactory<T>, bufferFactory: BufferFactory<T>,
vararg shape: Int, vararg shape: Int,
): BufferedNDSpace<T, A> = BufferedNDSpace(shape, space, bufferFactory) ): BufferedNDGroup<T, A> = BufferedNDGroup(shape, space, bufferFactory)
public inline fun <T, A : Space<T>, R> A.ndSpace( public inline fun <T, A : Group<T>, R> A.ndGroup(
noinline bufferFactory: BufferFactory<T>, noinline bufferFactory: BufferFactory<T>,
vararg shape: Int, vararg shape: Int,
action: BufferedNDSpace<T, A>.() -> R, action: BufferedNDGroup<T, A>.() -> R,
): R { ): R {
contract { callsInPlace(action, InvocationKind.EXACTLY_ONCE) } contract { callsInPlace(action, InvocationKind.EXACTLY_ONCE) }
return NDAlgebra.space(this, bufferFactory, *shape).run(action) return NDAlgebra.group(this, bufferFactory, *shape).run(action)
} }
//ring factories //ring factories

@ -1,8 +1,6 @@
package space.kscience.kmath.nd package space.kscience.kmath.nd
import space.kscience.kmath.operations.Field import space.kscience.kmath.operations.*
import space.kscience.kmath.operations.Ring
import space.kscience.kmath.operations.Space
import space.kscience.kmath.structures.* import space.kscience.kmath.structures.*
/** /**
@ -21,7 +19,7 @@ public class ShapeMismatchException(public val expected: IntArray, public val ac
* @param C the type of the element context. * @param C the type of the element context.
* @param N the type of the structure. * @param N the type of the structure.
*/ */
public interface NDAlgebra<T, C> { public interface NDAlgebra<T, C: Algebra<T>> {
/** /**
* The shape of ND-structures this algebra operates on. * The shape of ND-structures this algebra operates on.
*/ */
@ -33,7 +31,7 @@ public interface NDAlgebra<T, C> {
public val elementContext: C public val elementContext: C
/** /**
* Produces a new [N] structure using given initializer function. * Produces a new NDStructure using given initializer function.
*/ */
public fun produce(initializer: C.(IntArray) -> T): NDStructure<T> public fun produce(initializer: C.(IntArray) -> T): NDStructure<T>
@ -56,7 +54,7 @@ public interface NDAlgebra<T, C> {
* Element-wise invocation of function working on [T] on a [NDStructure]. * Element-wise invocation of function working on [T] on a [NDStructure].
*/ */
public operator fun Function1<T, T>.invoke(structure: NDStructure<T>): NDStructure<T> = public operator fun Function1<T, T>.invoke(structure: NDStructure<T>): NDStructure<T> =
structure.map() { value -> this@invoke(value) } structure.map { value -> this@invoke(value) }
public companion object public companion object
} }
@ -67,7 +65,7 @@ public interface NDAlgebra<T, C> {
* @param structures the structures to check. * @param structures the structures to check.
* @return the array of valid structures. * @return the array of valid structures.
*/ */
internal fun <T, C> NDAlgebra<T, C>.checkShape(vararg structures: NDStructure<T>): Array<out NDStructure<T>> = structures internal fun <T, C: Algebra<T>> NDAlgebra<T, C>.checkShape(vararg structures: NDStructure<T>): Array<out NDStructure<T>> = structures
.map(NDStructure<T>::shape) .map(NDStructure<T>::shape)
.singleOrNull { !shape.contentEquals(it) } .singleOrNull { !shape.contentEquals(it) }
?.let<IntArray, Array<out NDStructure<T>>> { throw ShapeMismatchException(shape, it) } ?.let<IntArray, Array<out NDStructure<T>>> { throw ShapeMismatchException(shape, it) }
@ -79,7 +77,7 @@ internal fun <T, C> NDAlgebra<T, C>.checkShape(vararg structures: NDStructure<T>
* @param element the structure to check. * @param element the structure to check.
* @return the valid structure. * @return the valid structure.
*/ */
internal fun <T, C> NDAlgebra<T, C>.checkShape(element: NDStructure<T>): NDStructure<T> { internal fun <T, C: Algebra<T>> NDAlgebra<T, C>.checkShape(element: NDStructure<T>): NDStructure<T> {
if (!element.shape.contentEquals(shape)) throw ShapeMismatchException(shape, element.shape) if (!element.shape.contentEquals(shape)) throw ShapeMismatchException(shape, element.shape)
return element return element
} }
@ -91,7 +89,7 @@ internal fun <T, C> NDAlgebra<T, C>.checkShape(element: NDStructure<T>): NDStruc
* @param N the type of ND structure. * @param N the type of ND structure.
* @param S the type of space of structure elements. * @param S the type of space of structure elements.
*/ */
public interface NDSpace<T, S : Space<T>> : Space<NDStructure<T>>, NDAlgebra<T, S> { public interface NDGroup<T, S : Group<T>> : Group<NDStructure<T>>, NDAlgebra<T, S> {
/** /**
* Element-wise addition. * Element-wise addition.
* *
@ -102,14 +100,14 @@ public interface NDSpace<T, S : Space<T>> : Space<NDStructure<T>>, NDAlgebra<T,
public override fun add(a: NDStructure<T>, b: NDStructure<T>): NDStructure<T> = public override fun add(a: NDStructure<T>, b: NDStructure<T>): NDStructure<T> =
combine(a, b) { aValue, bValue -> add(aValue, bValue) } combine(a, b) { aValue, bValue -> add(aValue, bValue) }
/** // /**
* Element-wise multiplication by scalar. // * Element-wise multiplication by scalar.
* // *
* @param a the multiplicand. // * @param a the multiplicand.
* @param k the multiplier. // * @param k the multiplier.
* @return the product. // * @return the product.
*/ // */
public override fun multiply(a: NDStructure<T>, k: Number): NDStructure<T> = a.map() { multiply(it, k) } // public override fun multiply(a: NDStructure<T>, k: Number): NDStructure<T> = a.map { multiply(it, k) }
// TODO move to extensions after KEEP-176 // TODO move to extensions after KEEP-176
@ -120,7 +118,7 @@ public interface NDSpace<T, S : Space<T>> : Space<NDStructure<T>>, NDAlgebra<T,
* @param arg the augend. * @param arg the augend.
* @return the sum. * @return the sum.
*/ */
public operator fun NDStructure<T>.plus(arg: T): NDStructure<T> = this.map() { value -> add(arg, value) } public operator fun NDStructure<T>.plus(arg: T): NDStructure<T> = this.map { value -> add(arg, value) }
/** /**
* Subtracts an element from ND structure of it. * Subtracts an element from ND structure of it.
@ -129,7 +127,7 @@ public interface NDSpace<T, S : Space<T>> : Space<NDStructure<T>>, NDAlgebra<T,
* @param arg the divisor. * @param arg the divisor.
* @return the quotient. * @return the quotient.
*/ */
public operator fun NDStructure<T>.minus(arg: T): NDStructure<T> = this.map() { value -> add(arg, -value) } public operator fun NDStructure<T>.minus(arg: T): NDStructure<T> = this.map { value -> add(arg, -value) }
/** /**
* Adds an element to ND structure of it. * Adds an element to ND structure of it.
@ -138,7 +136,7 @@ public interface NDSpace<T, S : Space<T>> : Space<NDStructure<T>>, NDAlgebra<T,
* @param arg the augend. * @param arg the augend.
* @return the sum. * @return the sum.
*/ */
public operator fun T.plus(arg: NDStructure<T>): NDStructure<T> = arg.map() { value -> add(this@plus, value) } public operator fun T.plus(arg: NDStructure<T>): NDStructure<T> = arg.map { value -> add(this@plus, value) }
/** /**
* Subtracts an ND structure from an element of it. * Subtracts an ND structure from an element of it.
@ -147,7 +145,7 @@ public interface NDSpace<T, S : Space<T>> : Space<NDStructure<T>>, NDAlgebra<T,
* @param arg the divisor. * @param arg the divisor.
* @return the quotient. * @return the quotient.
*/ */
public operator fun T.minus(arg: NDStructure<T>): NDStructure<T> = arg.map() { value -> add(-this@minus, value) } public operator fun T.minus(arg: NDStructure<T>): NDStructure<T> = arg.map { value -> add(-this@minus, value) }
public companion object public companion object
} }
@ -159,7 +157,7 @@ public interface NDSpace<T, S : Space<T>> : Space<NDStructure<T>>, NDAlgebra<T,
* @param N the type of ND structure. * @param N the type of ND structure.
* @param R the type of ring of structure elements. * @param R the type of ring of structure elements.
*/ */
public interface NDRing<T, R : Ring<T>> : Ring<NDStructure<T>>, NDSpace<T, R> { public interface NDRing<T, R : Ring<T>> : Ring<NDStructure<T>>, NDGroup<T, R> {
/** /**
* Element-wise multiplication. * Element-wise multiplication.
* *
@ -179,7 +177,7 @@ public interface NDRing<T, R : Ring<T>> : Ring<NDStructure<T>>, NDSpace<T, R> {
* @param arg the multiplier. * @param arg the multiplier.
* @return the product. * @return the product.
*/ */
public operator fun NDStructure<T>.times(arg: T): NDStructure<T> = this.map() { value -> multiply(arg, value) } public operator fun NDStructure<T>.times(arg: T): NDStructure<T> = this.map { value -> multiply(arg, value) }
/** /**
* Multiplies an element by a ND structure of it. * Multiplies an element by a ND structure of it.
@ -188,7 +186,7 @@ public interface NDRing<T, R : Ring<T>> : Ring<NDStructure<T>>, NDSpace<T, R> {
* @param arg the multiplier. * @param arg the multiplier.
* @return the product. * @return the product.
*/ */
public operator fun T.times(arg: NDStructure<T>): NDStructure<T> = arg.map() { value -> multiply(this@times, value) } public operator fun T.times(arg: NDStructure<T>): NDStructure<T> = arg.map { value -> multiply(this@times, value) }
public companion object public companion object
} }
@ -200,7 +198,7 @@ public interface NDRing<T, R : Ring<T>> : Ring<NDStructure<T>>, NDSpace<T, R> {
* @param N the type of ND structure. * @param N the type of ND structure.
* @param F the type field of structure elements. * @param F the type field of structure elements.
*/ */
public interface NDField<T, F : Field<T>> : Field<NDStructure<T>>, NDRing<T, F> { public interface NDField<T, F : Field<T>> : Field<NDStructure<T>>, NDRing<T, F>, ScaleOperations<NDStructure<T>> {
/** /**
* Element-wise division. * Element-wise division.
* *
@ -219,7 +217,7 @@ public interface NDField<T, F : Field<T>> : Field<NDStructure<T>>, NDRing<T, F>
* @param arg the divisor. * @param arg the divisor.
* @return the quotient. * @return the quotient.
*/ */
public operator fun NDStructure<T>.div(arg: T): NDStructure<T> = this.map() { value -> divide(arg, value) } public operator fun NDStructure<T>.div(arg: T): NDStructure<T> = this.map { value -> divide(arg, value) }
/** /**
* Divides an element by an ND structure of it. * Divides an element by an ND structure of it.
@ -228,7 +226,7 @@ public interface NDField<T, F : Field<T>> : Field<NDStructure<T>>, NDRing<T, F>
* @param arg the divisor. * @param arg the divisor.
* @return the quotient. * @return the quotient.
*/ */
public operator fun T.div(arg: NDStructure<T>): NDStructure<T> = arg.map() { divide(it, this@div) } public operator fun T.div(arg: NDStructure<T>): NDStructure<T> = arg.map { divide(it, this@div) }
// @ThreadLocal // @ThreadLocal
// public companion object { // public companion object {

@ -52,7 +52,7 @@ public interface NDStructure<T> {
* optimize operations and performance. If the feature is not present, null is defined. * optimize operations and performance. If the feature is not present, null is defined.
*/ */
@UnstableKMathAPI @UnstableKMathAPI
public fun <T : Any> getFeature(type: KClass<T>): T? = null public fun <F : Any> getFeature(type: KClass<F>): F? = null
public companion object { public companion object {
/** /**
@ -74,7 +74,7 @@ public interface NDStructure<T> {
* *
* Strides should be reused if possible. * Strides should be reused if possible.
*/ */
public fun <T> build( public fun <T> buffered(
strides: Strides, strides: Strides,
bufferFactory: BufferFactory<T> = Buffer.Companion::boxing, bufferFactory: BufferFactory<T> = Buffer.Companion::boxing,
initializer: (IntArray) -> T, initializer: (IntArray) -> T,
@ -94,11 +94,11 @@ public interface NDStructure<T> {
crossinline initializer: (IntArray) -> T, crossinline initializer: (IntArray) -> T,
): NDBuffer<T> = NDBuffer(strides, Buffer.auto(type, strides.linearSize) { i -> initializer(strides.index(i)) }) ): NDBuffer<T> = NDBuffer(strides, Buffer.auto(type, strides.linearSize) { i -> initializer(strides.index(i)) })
public fun <T> build( public fun <T> buffered(
shape: IntArray, shape: IntArray,
bufferFactory: BufferFactory<T> = Buffer.Companion::boxing, bufferFactory: BufferFactory<T> = Buffer.Companion::boxing,
initializer: (IntArray) -> T, initializer: (IntArray) -> T,
): NDBuffer<T> = build(DefaultStrides(shape), bufferFactory, initializer) ): NDBuffer<T> = buffered(DefaultStrides(shape), bufferFactory, initializer)
public inline fun <reified T : Any> auto( public inline fun <reified T : Any> auto(
shape: IntArray, shape: IntArray,

@ -2,8 +2,9 @@ package space.kscience.kmath.nd
import space.kscience.kmath.misc.UnstableKMathAPI import space.kscience.kmath.misc.UnstableKMathAPI
import space.kscience.kmath.operations.ExtendedField import space.kscience.kmath.operations.ExtendedField
import space.kscience.kmath.operations.NumbersAddOperations
import space.kscience.kmath.operations.RealField import space.kscience.kmath.operations.RealField
import space.kscience.kmath.operations.RingWithNumbers import space.kscience.kmath.operations.ScaleOperations
import space.kscience.kmath.structures.Buffer import space.kscience.kmath.structures.Buffer
import space.kscience.kmath.structures.RealBuffer import space.kscience.kmath.structures.RealBuffer
import kotlin.contracts.InvocationKind import kotlin.contracts.InvocationKind
@ -13,7 +14,8 @@ import kotlin.contracts.contract
public class RealNDField( public class RealNDField(
shape: IntArray, shape: IntArray,
) : BufferedNDField<Double, RealField>(shape, RealField, Buffer.Companion::real), ) : BufferedNDField<Double, RealField>(shape, RealField, Buffer.Companion::real),
RingWithNumbers<NDStructure<Double>>, NumbersAddOperations<NDStructure<Double>>,
ScaleOperations<NDStructure<Double>>,
ExtendedField<NDStructure<Double>> { ExtendedField<NDStructure<Double>> {
override val zero: NDBuffer<Double> by lazy { produce { zero } } override val zero: NDBuffer<Double> by lazy { produce { zero } }
@ -75,6 +77,8 @@ public class RealNDField(
return NDBuffer(strides, buffer) return NDBuffer(strides, buffer)
} }
override fun scale(a: NDStructure<Double>, value: Double): NDStructure<Double> = a.map { it * value }
override fun power(arg: NDStructure<Double>, pow: Number): NDBuffer<Double> = arg.map { power(it, pow) } override fun power(arg: NDStructure<Double>, pow: Number): NDBuffer<Double> = arg.map { power(it, pow) }
override fun exp(arg: NDStructure<Double>): NDBuffer<Double> = arg.map { exp(it) } override fun exp(arg: NDStructure<Double>): NDBuffer<Double> = arg.map { exp(it) }

@ -1,7 +1,7 @@
package space.kscience.kmath.nd package space.kscience.kmath.nd
import space.kscience.kmath.misc.UnstableKMathAPI import space.kscience.kmath.misc.UnstableKMathAPI
import space.kscience.kmath.operations.RingWithNumbers import space.kscience.kmath.operations.NumbersAddOperations
import space.kscience.kmath.operations.ShortRing import space.kscience.kmath.operations.ShortRing
import space.kscience.kmath.structures.Buffer import space.kscience.kmath.structures.Buffer
import space.kscience.kmath.structures.ShortBuffer import space.kscience.kmath.structures.ShortBuffer
@ -12,7 +12,7 @@ import kotlin.contracts.contract
public class ShortNDRing( public class ShortNDRing(
shape: IntArray, shape: IntArray,
) : BufferedNDRing<Short, ShortRing>(shape, ShortRing, Buffer.Companion::auto), ) : BufferedNDRing<Short, ShortRing>(shape, ShortRing, Buffer.Companion::auto),
RingWithNumbers<NDStructure<Short>> { NumbersAddOperations<NDStructure<Short>> {
override val zero: NDBuffer<Short> by lazy { produce { zero } } override val zero: NDBuffer<Short> by lazy { produce { zero } }
override val one: NDBuffer<Short> by lazy { produce { one } } override val one: NDBuffer<Short> by lazy { produce { one } }

@ -45,10 +45,12 @@ private inline class Buffer1DWrapper<T>(val buffer: Buffer<T>) : Structure1D<T>
/** /**
* Represent a [NDStructure] as [Structure1D]. Throw error in case of dimension mismatch * Represent a [NDStructure] as [Structure1D]. Throw error in case of dimension mismatch
*/ */
public fun <T> NDStructure<T>.as1D(): Structure1D<T> = if (shape.size == 1) { public fun <T> NDStructure<T>.as1D(): Structure1D<T> = this as? Structure1D<T> ?: if (shape.size == 1) {
if (this is NDBuffer) Buffer1DWrapper(this.buffer) else Structure1DWrapper(this) when (this) {
} else is NDBuffer -> Buffer1DWrapper(this.buffer)
error("Can't create 1d-structure from ${shape.size}d-structure") else -> Structure1DWrapper(this)
}
} else error("Can't create 1d-structure from ${shape.size}d-structure")
/** /**
* Represent this buffer as 1D structure * Represent this buffer as 1D structure

@ -1,9 +1,9 @@
package space.kscience.kmath.nd package space.kscience.kmath.nd
import space.kscience.kmath.linear.BufferMatrix import space.kscience.kmath.misc.UnstableKMathAPI
import space.kscience.kmath.linear.RealMatrixContext
import space.kscience.kmath.structures.Buffer import space.kscience.kmath.structures.Buffer
import space.kscience.kmath.structures.VirtualBuffer import space.kscience.kmath.structures.VirtualBuffer
import kotlin.reflect.KClass
/** /**
* A structure that is guaranteed to be two-dimensional. * A structure that is guaranteed to be two-dimensional.
@ -26,14 +26,14 @@ public interface Structure2D<T> : NDStructure<T> {
/** /**
* The buffer of rows of this structure. It gets elements from the structure dynamically. * The buffer of rows of this structure. It gets elements from the structure dynamically.
*/ */
public val rows: Buffer<Buffer<T>> public val rows: List<Buffer<T>>
get() = VirtualBuffer(rowNum) { i -> VirtualBuffer(colNum) { j -> get(i, j) } } get() = List(rowNum) { i -> VirtualBuffer(colNum) { j -> get(i, j) } }
/** /**
* The buffer of columns of this structure. It gets elements from the structure dynamically. * The buffer of columns of this structure. It gets elements from the structure dynamically.
*/ */
public val columns: Buffer<Buffer<T>> public val columns: List<Buffer<T>>
get() = VirtualBuffer(colNum) { j -> VirtualBuffer(rowNum) { i -> get(i, j) } } get() = List(colNum) { j -> VirtualBuffer(rowNum) { i -> get(i, j) } }
/** /**
* Retrieves an element from the structure by two indices. * Retrieves an element from the structure by two indices.
@ -54,21 +54,13 @@ public interface Structure2D<T> : NDStructure<T> {
for (j in 0 until colNum) yield(intArrayOf(i, j) to get(i, j)) for (j in 0 until colNum) yield(intArrayOf(i, j) to get(i, j))
} }
public companion object { public companion object
public inline fun real(
rows: Int,
columns: Int,
crossinline init: (i: Int, j: Int) -> Double,
): BufferMatrix<Double> = RealMatrixContext.produce(rows,columns) { i, j ->
init(i, j)
}
}
} }
/** /**
* A 2D wrapper for nd-structure * A 2D wrapper for nd-structure
*/ */
private inline class Structure2DWrapper<T>(val structure: NDStructure<T>) : Structure2D<T> { private class Structure2DWrapper<T>(val structure: NDStructure<T>) : Structure2D<T> {
override val shape: IntArray get() = structure.shape override val shape: IntArray get() = structure.shape
override val rowNum: Int get() = shape[0] override val rowNum: Int get() = shape[0]
@ -76,20 +68,22 @@ private inline class Structure2DWrapper<T>(val structure: NDStructure<T>) : Stru
override operator fun get(i: Int, j: Int): T = structure[i, j] override operator fun get(i: Int, j: Int): T = structure[i, j]
@UnstableKMathAPI
override fun <F : Any> getFeature(type: KClass<F>): F? = structure.getFeature(type)
override fun elements(): Sequence<Pair<IntArray, T>> = structure.elements() override fun elements(): Sequence<Pair<IntArray, T>> = structure.elements()
override fun equals(other: Any?): Boolean = structure == other
override fun hashCode(): Int = structure.hashCode()
} }
/** /**
* Represent a [NDStructure] as [Structure1D]. Throw error in case of dimension mismatch * Represent a [NDStructure] as [Structure1D]. Throw error in case of dimension mismatch
*/ */
public fun <T> NDStructure<T>.as2D(): Structure2D<T> = if (shape.size == 2) public fun <T> NDStructure<T>.as2D(): Structure2D<T> = this as? Structure2D<T> ?: when (shape.size) {
Structure2DWrapper(this) 2 -> Structure2DWrapper(this)
else else -> error("Can't create 2d-structure from ${shape.size}d-structure")
error("Can't create 2d-structure from ${shape.size}d-structure") }
/** internal fun <T> Structure2D<T>.unwrap(): NDStructure<T> = if (this is Structure2DWrapper) structure else this
* Alias for [Structure2D] with more familiar name.
*
* @param T the type of items in the matrix.
*/
public typealias Matrix<T> = Structure2D<T>

@ -87,10 +87,11 @@ public interface Algebra<T> {
* @param right the second argument of operation. * @param right the second argument of operation.
* @return a result of operation. * @return a result of operation.
*/ */
public fun binaryOperation(operation: String, left: T, right: T): T = binaryOperationFunction(operation)(left, right) public fun binaryOperation(operation: String, left: T, right: T): T =
binaryOperationFunction(operation)(left, right)
} }
public fun <T: Any> Algebra<T>.bindSymbol(symbol: Symbol): T = bindSymbol(symbol.identity) public fun <T : Any> Algebra<T>.bindSymbol(symbol: Symbol): T = bindSymbol(symbol.identity)
/** /**
* Call a block with an [Algebra] as receiver. * Call a block with an [Algebra] as receiver.
@ -104,7 +105,7 @@ public inline operator fun <A : Algebra<*>, R> A.invoke(block: A.() -> R): R = r
* *
* @param T the type of element of this semispace. * @param T the type of element of this semispace.
*/ */
public interface SpaceOperations<T> : Algebra<T> { public interface GroupOperations<T> : Algebra<T> {
/** /**
* Addition of two elements. * Addition of two elements.
* *
@ -114,15 +115,6 @@ public interface SpaceOperations<T> : Algebra<T> {
*/ */
public fun add(a: T, b: T): T public fun add(a: T, b: T): T
/**
* Multiplication of element by scalar.
*
* @param a the multiplier.
* @param k the multiplicand.
* @return the produce.
*/
public fun multiply(a: T, k: Number): T
// Operations to be performed in this context. Could be moved to extensions in case of KEEP-176 // Operations to be performed in this context. Could be moved to extensions in case of KEEP-176
/** /**
@ -131,7 +123,7 @@ public interface SpaceOperations<T> : Algebra<T> {
* @receiver this value. * @receiver this value.
* @return the additive inverse of this value. * @return the additive inverse of this value.
*/ */
public operator fun T.unaryMinus(): T = multiply(this, -1.0) public operator fun T.unaryMinus(): T
/** /**
* Returns this value. * Returns this value.
@ -159,36 +151,8 @@ public interface SpaceOperations<T> : Algebra<T> {
*/ */
public operator fun T.minus(b: T): T = add(this, -b) public operator fun T.minus(b: T): T = add(this, -b)
/**
* Multiplication of this element by a scalar.
*
* @receiver the multiplier.
* @param k the multiplicand.
* @return the product.
*/
public operator fun T.times(k: Number): T = multiply(this, k)
/**
* Division of this element by scalar.
*
* @receiver the dividend.
* @param k the divisor.
* @return the quotient.
*/
@Deprecated("Dividing not allowed in a Ring")
public operator fun T.div(k: Number): T = multiply(this, 1.0 / k.toDouble())
/**
* Multiplication of this number by element.
*
* @receiver the multiplier.
* @param b the multiplicand.
* @return the product.
*/
public operator fun Number.times(b: T): T = b * this
public override fun unaryOperationFunction(operation: String): (arg: T) -> T = when (operation) { public override fun unaryOperationFunction(operation: String): (arg: T) -> T = when (operation) {
PLUS_OPERATION -> { arg -> arg } PLUS_OPERATION -> { arg -> +arg }
MINUS_OPERATION -> { arg -> -arg } MINUS_OPERATION -> { arg -> -arg }
else -> super.unaryOperationFunction(operation) else -> super.unaryOperationFunction(operation)
} }
@ -213,12 +177,11 @@ public interface SpaceOperations<T> : Algebra<T> {
} }
/** /**
* Represents linear space with neutral element, i.e. algebraic structure with associative, binary operation [add] and * Represents linear space with neutral element, i.e. algebraic structure with associative, binary operation [add].
* scalar multiplication [multiply].
* *
* @param T the type of element of this semispace. * @param T the type of element of this semispace.
*/ */
public interface Space<T> : SpaceOperations<T> { public interface Group<T> : GroupOperations<T> {
/** /**
* The neutral element of addition. * The neutral element of addition.
*/ */
@ -231,7 +194,7 @@ public interface Space<T> : SpaceOperations<T> {
* *
* @param T the type of element of this semiring. * @param T the type of element of this semiring.
*/ */
public interface RingOperations<T> : SpaceOperations<T> { public interface RingOperations<T> : GroupOperations<T> {
/** /**
* Multiplies two elements. * Multiplies two elements.
* *
@ -267,7 +230,7 @@ public interface RingOperations<T> : SpaceOperations<T> {
* *
* @param T the type of element of this ring. * @param T the type of element of this ring.
*/ */
public interface Ring<T> : Space<T>, RingOperations<T> { public interface Ring<T> : Group<T>, RingOperations<T> {
/** /**
* neutral operation for multiplication * neutral operation for multiplication
*/ */
@ -318,13 +281,6 @@ public interface FieldOperations<T> : RingOperations<T> {
* *
* @param T the type of element of this semifield. * @param T the type of element of this semifield.
*/ */
public interface Field<T> : Ring<T>, FieldOperations<T> { public interface Field<T> : Ring<T>, FieldOperations<T>, ScaleOperations<T>, NumericAlgebra<T> {
/** override fun number(value: Number): T = scale(one, value.toDouble())
* Division of element by scalar. }
*
* @receiver the dividend.
* @param b the divisor.
* @return the quotient.
*/
public operator fun Number.div(b: T): T = this * divide(one, b)
}

@ -14,24 +14,24 @@ public interface AlgebraElement<T, C : Algebra<T>> {
*/ */
public val context: C public val context: C
} }
//
/** ///**
* Divides this element by number. // * Divides this element by number.
* // *
* @param k the divisor. // * @param k the divisor.
* @return the quotient. // * @return the quotient.
*/ // */
public operator fun <T : AlgebraElement<T, S>, S : Space<T>> T.div(k: Number): T = //public operator fun <T : AlgebraElement<T, S>, S : Space<T>> T.div(k: Number): T =
context.multiply(this, 1.0 / k.toDouble()) // context.multiply(this, 1.0 / k.toDouble())
//
/** ///**
* Multiplies this element by number. // * Multiplies this element by number.
* // *
* @param k the multiplicand. // * @param k the multiplicand.
* @return the product. // * @return the product.
*/ // */
public operator fun <T : AlgebraElement<T, S>, S : Space<T>> T.times(k: Number): T = //public operator fun <T : AlgebraElement<T, S>, S : Space<T>> T.times(k: Number): T =
context.multiply(this, k.toDouble()) // context.multiply(this, k.toDouble())
/** /**
* Subtracts element from this one. * Subtracts element from this one.
@ -39,8 +39,9 @@ public operator fun <T : AlgebraElement<T, S>, S : Space<T>> T.times(k: Number):
* @param b the subtrahend. * @param b the subtrahend.
* @return the difference. * @return the difference.
*/ */
public operator fun <T : AlgebraElement<T, S>, S : Space<T>> T.minus(b: T): T = @UnstableKMathAPI
context.add(this, context.multiply(b, -1.0)) public operator fun <T : AlgebraElement<T, S>, S : NumbersAddOperations<T>> T.minus(b: T): T =
context.add(this, context.run { -b})
/** /**
* Adds element to this one. * Adds element to this one.
@ -48,14 +49,14 @@ public operator fun <T : AlgebraElement<T, S>, S : Space<T>> T.minus(b: T): T =
* @param b the augend. * @param b the augend.
* @return the sum. * @return the sum.
*/ */
public operator fun <T : AlgebraElement<T, S>, S : Space<T>> T.plus(b: T): T = public operator fun <T : AlgebraElement<T, S>, S : Group<T>> T.plus(b: T): T =
context.add(this, b) context.add(this, b)
/** ///**
* Number times element // * Number times element
*/ // */
public operator fun <T : AlgebraElement<T, S>, S : Space<T>> Number.times(element: T): T = //public operator fun <T : AlgebraElement<T, S>, S : Space<T>> Number.times(element: T): T =
element.times(this) // element.times(this)
/** /**
@ -79,14 +80,14 @@ public operator fun <T : AlgebraElement<T, F>, F : Field<T>> T.div(b: T): T =
/** /**
* The element of [Space]. * The element of [Group].
* *
* @param T the type of space operation results. * @param T the type of space operation results.
* @param I self type of the element. Needed for static type checking. * @param I self type of the element. Needed for static type checking.
* @param S the type of space. * @param S the type of space.
*/ */
@UnstableKMathAPI @UnstableKMathAPI
public interface SpaceElement<T : SpaceElement<T, S>, S : Space<T>> : AlgebraElement<T, S> public interface SpaceElement<T : SpaceElement<T, S>, S : Group<T>> : AlgebraElement<T, S>
/** /**
* The element of [Ring]. * The element of [Ring].

@ -21,29 +21,28 @@ public typealias TBase = ULong
* @author Robert Drynkin (https://github.com/robdrynkin) and Peter Klimai (https://github.com/pklimai) * @author Robert Drynkin (https://github.com/robdrynkin) and Peter Klimai (https://github.com/pklimai)
*/ */
@OptIn(UnstableKMathAPI::class) @OptIn(UnstableKMathAPI::class)
public object BigIntField : Field<BigInt>, RingWithNumbers<BigInt> { public object BigIntField : Field<BigInt>, NumbersAddOperations<BigInt>, ScaleOperations<BigInt> {
override val zero: BigInt = BigInt.ZERO override val zero: BigInt = BigInt.ZERO
override val one: BigInt = BigInt.ONE override val one: BigInt = BigInt.ONE
override fun add(a: BigInt, b: BigInt): BigInt = a.plus(b)
override fun number(value: Number): BigInt = value.toLong().toBigInt() override fun number(value: Number): BigInt = value.toLong().toBigInt()
override fun multiply(a: BigInt, k: Number): BigInt = a.times(number(k)) @Suppress("EXTENSION_SHADOWED_BY_MEMBER")
override fun BigInt.unaryMinus(): BigInt = -this
override fun add(a: BigInt, b: BigInt): BigInt = a.plus(b)
override fun scale(a: BigInt, value: Double): BigInt = a.times(number(value))
override fun multiply(a: BigInt, b: BigInt): BigInt = a.times(b) override fun multiply(a: BigInt, b: BigInt): BigInt = a.times(b)
override fun divide(a: BigInt, b: BigInt): BigInt = a.div(b)
public operator fun String.unaryPlus(): BigInt = this.parseBigInteger() ?: error("Can't parse $this as big integer") public operator fun String.unaryPlus(): BigInt = this.parseBigInteger() ?: error("Can't parse $this as big integer")
public operator fun String.unaryMinus(): BigInt = public operator fun String.unaryMinus(): BigInt =
-(this.parseBigInteger() ?: error("Can't parse $this as big integer")) -(this.parseBigInteger() ?: error("Can't parse $this as big integer"))
override fun divide(a: BigInt, b: BigInt): BigInt = a.div(b)
} }
public class BigInt internal constructor( public class BigInt internal constructor(
private val sign: Byte, private val sign: Byte,
private val magnitude: Magnitude private val magnitude: Magnitude,
) : Comparable<BigInt> { ) : Comparable<BigInt> {
public override fun compareTo(other: BigInt): Int = when { public override fun compareTo(other: BigInt): Int = when {
(sign == 0.toByte()) and (other.sign == 0.toByte()) -> 0 (sign == 0.toByte()) and (other.sign == 0.toByte()) -> 0
sign < other.sign -> -1 sign < other.sign -> -1

@ -81,14 +81,53 @@ public interface NumericAlgebra<T> : Algebra<T> {
rightSideNumberOperationFunction(operation)(left, right) rightSideNumberOperationFunction(operation)(left, right)
} }
/**
* Scale by scalar operations
*/
public interface ScaleOperations<T> : Algebra<T> {
/**
* Scaling an element by a scalar.
*
* @param a the multiplier.
* @param value the multiplicand.
* @return the produce.
*/
public fun scale(a: T, value: Double): T
/**
* Multiplication of this element by a scalar.
*
* @receiver the multiplier.
* @param k the multiplicand.
* @return the product.
*/
public operator fun T.times(k: Number): T = scale(this, k.toDouble())
/**
* Division of this element by scalar.
*
* @receiver the dividend.
* @param k the divisor.
* @return the quotient.
*/
public operator fun T.div(k: Number): T = scale(this, 1.0 / k.toDouble())
/**
* Multiplication of this number by element.
*
* @receiver the multiplier.
* @param b the multiplicand.
* @return the product.
*/
public operator fun Number.times(b: T): T = b * this
}
/** /**
* A combination of [NumericAlgebra] and [Ring] that adds intrinsic simple operations on numbers like `T+1` * A combination of [NumericAlgebra] and [Ring] that adds intrinsic simple operations on numbers like `T+1`
* TODO to be removed and replaced by extensions after multiple receivers are there * TODO to be removed and replaced by extensions after multiple receivers are there
*/ */
@UnstableKMathAPI @UnstableKMathAPI
public interface RingWithNumbers<T>: Ring<T>, NumericAlgebra<T>{ public interface NumbersAddOperations<T> : Group<T>, NumericAlgebra<T> {
public override fun number(value: Number): T = one * value
/** /**
* Addition of element and scalar. * Addition of element and scalar.
* *

@ -107,111 +107,6 @@ public fun <T : AlgebraElement<T, out TrigonometricOperations<T>>> acos(arg: T):
@UnstableKMathAPI @UnstableKMathAPI
public fun <T : AlgebraElement<T, out TrigonometricOperations<T>>> atan(arg: T): T = arg.context.atan(arg) public fun <T : AlgebraElement<T, out TrigonometricOperations<T>>> atan(arg: T): T = arg.context.atan(arg)
/**
* A container for hyperbolic trigonometric operations for specific type.
*
* @param T the type of element of this structure.
*/
public interface HyperbolicOperations<T> : Algebra<T> {
/**
* Computes the hyperbolic sine of [arg].
*/
public fun sinh(arg: T): T
/**
* Computes the hyperbolic cosine of [arg].
*/
public fun cosh(arg: T): T
/**
* Computes the hyperbolic tangent of [arg].
*/
public fun tanh(arg: T): T
/**
* Computes the inverse hyperbolic sine of [arg].
*/
public fun asinh(arg: T): T
/**
* Computes the inverse hyperbolic cosine of [arg].
*/
public fun acosh(arg: T): T
/**
* Computes the inverse hyperbolic tangent of [arg].
*/
public fun atanh(arg: T): T
public companion object {
/**
* The identifier of hyperbolic sine.
*/
public const val SINH_OPERATION: String = "sinh"
/**
* The identifier of hyperbolic cosine.
*/
public const val COSH_OPERATION: String = "cosh"
/**
* The identifier of hyperbolic tangent.
*/
public const val TANH_OPERATION: String = "tanh"
/**
* The identifier of inverse hyperbolic sine.
*/
public const val ASINH_OPERATION: String = "asinh"
/**
* The identifier of inverse hyperbolic cosine.
*/
public const val ACOSH_OPERATION: String = "acosh"
/**
* The identifier of inverse hyperbolic tangent.
*/
public const val ATANH_OPERATION: String = "atanh"
}
}
/**
* Computes the hyperbolic sine of [arg].
*/
@UnstableKMathAPI
public fun <T : AlgebraElement<T, out HyperbolicOperations<T>>> sinh(arg: T): T = arg.context.sinh(arg)
/**
* Computes the hyperbolic cosine of [arg].
*/
@UnstableKMathAPI
public fun <T : AlgebraElement<T, out HyperbolicOperations<T>>> cosh(arg: T): T = arg.context.cosh(arg)
/**
* Computes the hyperbolic tangent of [arg].
*/
@UnstableKMathAPI
public fun <T : AlgebraElement<T, out HyperbolicOperations<T>>> tanh(arg: T): T = arg.context.tanh(arg)
/**
* Computes the inverse hyperbolic sine of [arg].
*/
@UnstableKMathAPI
public fun <T : AlgebraElement<T, out HyperbolicOperations<T>>> asinh(arg: T): T = arg.context.asinh(arg)
/**
* Computes the inverse hyperbolic cosine of [arg].
*/
@UnstableKMathAPI
public fun <T : AlgebraElement<T, out HyperbolicOperations<T>>> acosh(arg: T): T = arg.context.acosh(arg)
/**
* Computes the inverse hyperbolic tangent of [arg].
*/
@UnstableKMathAPI
public fun <T : AlgebraElement<T, out HyperbolicOperations<T>>> atanh(arg: T): T = arg.context.atanh(arg)
/** /**
* A context extension to include power operations based on exponentiation. * A context extension to include power operations based on exponentiation.
* *
@ -284,6 +179,36 @@ public interface ExponentialOperations<T> : Algebra<T> {
*/ */
public fun ln(arg: T): T public fun ln(arg: T): T
/**
* Computes the hyperbolic sine of [arg].
*/
public fun sinh(arg: T): T
/**
* Computes the hyperbolic cosine of [arg].
*/
public fun cosh(arg: T): T
/**
* Computes the hyperbolic tangent of [arg].
*/
public fun tanh(arg: T): T
/**
* Computes the inverse hyperbolic sine of [arg].
*/
public fun asinh(arg: T): T
/**
* Computes the inverse hyperbolic cosine of [arg].
*/
public fun acosh(arg: T): T
/**
* Computes the inverse hyperbolic tangent of [arg].
*/
public fun atanh(arg: T): T
public companion object { public companion object {
/** /**
* The identifier of exponential function. * The identifier of exponential function.
@ -294,6 +219,36 @@ public interface ExponentialOperations<T> : Algebra<T> {
* The identifier of natural logarithm. * The identifier of natural logarithm.
*/ */
public const val LN_OPERATION: String = "ln" public const val LN_OPERATION: String = "ln"
/**
* The identifier of hyperbolic sine.
*/
public const val SINH_OPERATION: String = "sinh"
/**
* The identifier of hyperbolic cosine.
*/
public const val COSH_OPERATION: String = "cosh"
/**
* The identifier of hyperbolic tangent.
*/
public const val TANH_OPERATION: String = "tanh"
/**
* The identifier of inverse hyperbolic sine.
*/
public const val ASINH_OPERATION: String = "asinh"
/**
* The identifier of inverse hyperbolic cosine.
*/
public const val ACOSH_OPERATION: String = "acosh"
/**
* The identifier of inverse hyperbolic tangent.
*/
public const val ATANH_OPERATION: String = "atanh"
} }
} }
@ -309,6 +264,43 @@ public fun <T : AlgebraElement<T, out ExponentialOperations<T>>> exp(arg: T): T
@UnstableKMathAPI @UnstableKMathAPI
public fun <T : AlgebraElement<T, out ExponentialOperations<T>>> ln(arg: T): T = arg.context.ln(arg) public fun <T : AlgebraElement<T, out ExponentialOperations<T>>> ln(arg: T): T = arg.context.ln(arg)
/**
* Computes the hyperbolic sine of [arg].
*/
@UnstableKMathAPI
public fun <T : AlgebraElement<T, out ExponentialOperations<T>>> sinh(arg: T): T = arg.context.sinh(arg)
/**
* Computes the hyperbolic cosine of [arg].
*/
@UnstableKMathAPI
public fun <T : AlgebraElement<T, out ExponentialOperations<T>>> cosh(arg: T): T = arg.context.cosh(arg)
/**
* Computes the hyperbolic tangent of [arg].
*/
@UnstableKMathAPI
public fun <T : AlgebraElement<T, out ExponentialOperations<T>>> tanh(arg: T): T = arg.context.tanh(arg)
/**
* Computes the inverse hyperbolic sine of [arg].
*/
@UnstableKMathAPI
public fun <T : AlgebraElement<T, out ExponentialOperations<T>>> asinh(arg: T): T = arg.context.asinh(arg)
/**
* Computes the inverse hyperbolic cosine of [arg].
*/
@UnstableKMathAPI
public fun <T : AlgebraElement<T, out ExponentialOperations<T>>> acosh(arg: T): T = arg.context.acosh(arg)
/**
* Computes the inverse hyperbolic tangent of [arg].
*/
@UnstableKMathAPI
public fun <T : AlgebraElement<T, out ExponentialOperations<T>>> atanh(arg: T): T = arg.context.atanh(arg)
/** /**
* A container for norm functional on element. * A container for norm functional on element.
* *

@ -1,47 +1,49 @@
package space.kscience.kmath.operations package space.kscience.kmath.operations
/** /**
* Returns the sum of all elements in the iterable in this [Space]. * Returns the sum of all elements in the iterable in this [Group].
* *
* @receiver the algebra that provides addition. * @receiver the algebra that provides addition.
* @param data the iterable to sum up. * @param data the iterable to sum up.
* @return the sum. * @return the sum.
*/ */
public fun <T> Space<T>.sum(data: Iterable<T>): T = data.fold(zero) { left, right -> add(left, right) } public fun <T> Group<T>.sum(data: Iterable<T>): T = data.fold(zero) { left, right -> add(left, right) }
/** /**
* Returns the sum of all elements in the sequence in this [Space]. * Returns the sum of all elements in the sequence in this [Group].
* *
* @receiver the algebra that provides addition. * @receiver the algebra that provides addition.
* @param data the sequence to sum up. * @param data the sequence to sum up.
* @return the sum. * @return the sum.
*/ */
public fun <T> Space<T>.sum(data: Sequence<T>): T = data.fold(zero) { left, right -> add(left, right) } public fun <T> Group<T>.sum(data: Sequence<T>): T = data.fold(zero) { left, right -> add(left, right) }
/** /**
* Returns an average value of elements in the iterable in this [Space]. * Returns an average value of elements in the iterable in this [Group].
* *
* @receiver the algebra that provides addition and division. * @receiver the algebra that provides addition and division.
* @param data the iterable to find average. * @param data the iterable to find average.
* @return the average value. * @return the average value.
* @author Iaroslav Postovalov * @author Iaroslav Postovalov
*/ */
public fun <T> Space<T>.average(data: Iterable<T>): T = sum(data) / data.count() public fun <T, S> S.average(data: Iterable<T>): T where S : Group<T>, S : ScaleOperations<T> =
sum(data) / data.count()
/** /**
* Returns an average value of elements in the sequence in this [Space]. * Returns an average value of elements in the sequence in this [Group].
* *
* @receiver the algebra that provides addition and division. * @receiver the algebra that provides addition and division.
* @param data the sequence to find average. * @param data the sequence to find average.
* @return the average value. * @return the average value.
* @author Iaroslav Postovalov * @author Iaroslav Postovalov
*/ */
public fun <T> Space<T>.average(data: Sequence<T>): T = sum(data) / data.count() public fun <T, S> S.average(data: Sequence<T>): T where S : Group<T>, S : ScaleOperations<T> =
sum(data) / data.count()
/** /**
* Absolute of the comparable [value] * Absolute of the comparable [value]
*/ */
public fun <T : Comparable<T>> Space<T>.abs(value: T): T = if (value > zero) value else -value public fun <T : Comparable<T>> Group<T>.abs(value: T): T = if (value > zero) value else -value
/** /**
* Returns the sum of all elements in the iterable in provided space. * Returns the sum of all elements in the iterable in provided space.
@ -50,7 +52,7 @@ public fun <T : Comparable<T>> Space<T>.abs(value: T): T = if (value > zero) val
* @param space the algebra that provides addition. * @param space the algebra that provides addition.
* @return the sum. * @return the sum.
*/ */
public fun <T> Iterable<T>.sumWith(space: Space<T>): T = space.sum(this) public fun <T> Iterable<T>.sumWith(space: Group<T>): T = space.sum(this)
/** /**
* Returns the sum of all elements in the sequence in provided space. * Returns the sum of all elements in the sequence in provided space.
@ -59,27 +61,29 @@ public fun <T> Iterable<T>.sumWith(space: Space<T>): T = space.sum(this)
* @param space the algebra that provides addition. * @param space the algebra that provides addition.
* @return the sum. * @return the sum.
*/ */
public fun <T> Sequence<T>.sumWith(space: Space<T>): T = space.sum(this) public fun <T> Sequence<T>.sumWith(space: Group<T>): T = space.sum(this)
/** /**
* Returns an average value of elements in the iterable in this [Space]. * Returns an average value of elements in the iterable in this [Group].
* *
* @receiver the iterable to find average. * @receiver the iterable to find average.
* @param space the algebra that provides addition and division. * @param space the algebra that provides addition and division.
* @return the average value. * @return the average value.
* @author Iaroslav Postovalov * @author Iaroslav Postovalov
*/ */
public fun <T> Iterable<T>.averageWith(space: Space<T>): T = space.average(this) public fun <T, S> Iterable<T>.averageWith(space: S): T where S : Group<T>, S : ScaleOperations<T> =
space.average(this)
/** /**
* Returns an average value of elements in the sequence in this [Space]. * Returns an average value of elements in the sequence in this [Group].
* *
* @receiver the sequence to find average. * @receiver the sequence to find average.
* @param space the algebra that provides addition and division. * @param space the algebra that provides addition and division.
* @return the average value. * @return the average value.
* @author Iaroslav Postovalov * @author Iaroslav Postovalov
*/ */
public fun <T> Sequence<T>.averageWith(space: Space<T>): T = space.average(this) public fun <T, S> Sequence<T>.averageWith(space: S): T where S : Group<T>, S : ScaleOperations<T> =
space.average(this)
//TODO optimized power operation //TODO optimized power operation

@ -8,7 +8,6 @@ import kotlin.math.pow as kpow
public interface ExtendedFieldOperations<T> : public interface ExtendedFieldOperations<T> :
FieldOperations<T>, FieldOperations<T>,
TrigonometricOperations<T>, TrigonometricOperations<T>,
HyperbolicOperations<T>,
PowerOperations<T>, PowerOperations<T>,
ExponentialOperations<T> { ExponentialOperations<T> {
public override fun tan(arg: T): T = sin(arg) / cos(arg) public override fun tan(arg: T): T = sin(arg) / cos(arg)
@ -21,15 +20,15 @@ public interface ExtendedFieldOperations<T> :
TrigonometricOperations.ACOS_OPERATION -> ::acos TrigonometricOperations.ACOS_OPERATION -> ::acos
TrigonometricOperations.ASIN_OPERATION -> ::asin TrigonometricOperations.ASIN_OPERATION -> ::asin
TrigonometricOperations.ATAN_OPERATION -> ::atan TrigonometricOperations.ATAN_OPERATION -> ::atan
HyperbolicOperations.COSH_OPERATION -> ::cosh
HyperbolicOperations.SINH_OPERATION -> ::sinh
HyperbolicOperations.TANH_OPERATION -> ::tanh
HyperbolicOperations.ACOSH_OPERATION -> ::acosh
HyperbolicOperations.ASINH_OPERATION -> ::asinh
HyperbolicOperations.ATANH_OPERATION -> ::atanh
PowerOperations.SQRT_OPERATION -> ::sqrt PowerOperations.SQRT_OPERATION -> ::sqrt
ExponentialOperations.EXP_OPERATION -> ::exp ExponentialOperations.EXP_OPERATION -> ::exp
ExponentialOperations.LN_OPERATION -> ::ln ExponentialOperations.LN_OPERATION -> ::ln
ExponentialOperations.COSH_OPERATION -> ::cosh
ExponentialOperations.SINH_OPERATION -> ::sinh
ExponentialOperations.TANH_OPERATION -> ::tanh
ExponentialOperations.ACOSH_OPERATION -> ::acosh
ExponentialOperations.ASINH_OPERATION -> ::asinh
ExponentialOperations.ATANH_OPERATION -> ::atanh
else -> super<FieldOperations>.unaryOperationFunction(operation) else -> super<FieldOperations>.unaryOperationFunction(operation)
} }
} }
@ -37,18 +36,18 @@ public interface ExtendedFieldOperations<T> :
/** /**
* Advanced Number-like field that implements basic operations. * Advanced Number-like field that implements basic operations.
*/ */
public interface ExtendedField<T> : ExtendedFieldOperations<T>, Field<T>, NumericAlgebra<T> { public interface ExtendedField<T> : ExtendedFieldOperations<T>, Field<T>, NumericAlgebra<T>, ScaleOperations<T> {
public override fun sinh(arg: T): T = (exp(arg) - exp(-arg)) / 2 public override fun sinh(arg: T): T = (exp(arg) - exp(-arg)) / 2.0
public override fun cosh(arg: T): T = (exp(arg) + exp(-arg)) / 2 public override fun cosh(arg: T): T = (exp(arg) + exp(-arg)) / 2.0
public override fun tanh(arg: T): T = (exp(arg) - exp(-arg)) / (exp(-arg) + exp(arg)) public override fun tanh(arg: T): T = (exp(arg) - exp(-arg)) / (exp(-arg) + exp(arg))
public override fun asinh(arg: T): T = ln(sqrt(arg * arg + one) + arg) public override fun asinh(arg: T): T = ln(sqrt(arg * arg + one) + arg)
public override fun acosh(arg: T): T = ln(arg + sqrt((arg - one) * (arg + one))) public override fun acosh(arg: T): T = ln(arg + sqrt((arg - one) * (arg + one)))
public override fun atanh(arg: T): T = (ln(arg + one) - ln(one - arg)) / 2 public override fun atanh(arg: T): T = (ln(arg + one) - ln(one - arg)) / 2.0
public override fun rightSideNumberOperationFunction(operation: String): (left: T, right: Number) -> T = public override fun rightSideNumberOperationFunction(operation: String): (left: T, right: Number) -> T =
when (operation) { when (operation) {
PowerOperations.POW_OPERATION -> ::power PowerOperations.POW_OPERATION -> ::power
else -> super.rightSideNumberOperationFunction(operation) else -> super<Field>.rightSideNumberOperationFunction(operation)
} }
} }
@ -56,28 +55,27 @@ public interface ExtendedField<T> : ExtendedFieldOperations<T>, Field<T>, Numeri
* A field for [Double] without boxing. Does not produce appropriate field element. * A field for [Double] without boxing. Does not produce appropriate field element.
*/ */
@Suppress("EXTENSION_SHADOWED_BY_MEMBER", "OVERRIDE_BY_INLINE", "NOTHING_TO_INLINE") @Suppress("EXTENSION_SHADOWED_BY_MEMBER", "OVERRIDE_BY_INLINE", "NOTHING_TO_INLINE")
public object RealField : ExtendedField<Double>, Norm<Double, Double> { public object RealField : ExtendedField<Double>, Norm<Double, Double>, ScaleOperations<Double> {
public override val zero: Double public override val zero: Double = 0.0
get() = 0.0 public override val one: Double = 1.0
public override val one: Double
get() = 1.0
override fun number(value: Number): Double = value.toDouble() override fun number(value: Number): Double = value.toDouble()
public override fun binaryOperationFunction(operation: String): (left: Double, right: Double) -> Double = public override fun binaryOperationFunction(operation: String): (left: Double, right: Double) -> Double =
when (operation) { when (operation) {
PowerOperations.POW_OPERATION -> ::power PowerOperations.POW_OPERATION -> ::power
else -> super.binaryOperationFunction(operation) else -> super<ExtendedField>.binaryOperationFunction(operation)
} }
public override inline fun add(a: Double, b: Double): Double = a + b public override inline fun add(a: Double, b: Double): Double = a + b
public override inline fun multiply(a: Double, k: Number): Double = a * k.toDouble() // public override inline fun multiply(a: Double, k: Number): Double = a * k.toDouble()
// override fun divide(a: Double, k: Number): Double = a / k.toDouble()
public override inline fun multiply(a: Double, b: Double): Double = a * b public override inline fun multiply(a: Double, b: Double): Double = a * b
public override inline fun divide(a: Double, b: Double): Double = a / b public override inline fun divide(a: Double, b: Double): Double = a / b
override fun scale(a: Double, value: Double): Double = a * value
public override inline fun sin(arg: Double): Double = kotlin.math.sin(arg) public override inline fun sin(arg: Double): Double = kotlin.math.sin(arg)
public override inline fun cos(arg: Double): Double = kotlin.math.cos(arg) public override inline fun cos(arg: Double): Double = kotlin.math.cos(arg)
public override inline fun tan(arg: Double): Double = kotlin.math.tan(arg) public override inline fun tan(arg: Double): Double = kotlin.math.tan(arg)
@ -110,11 +108,8 @@ public object RealField : ExtendedField<Double>, Norm<Double, Double> {
*/ */
@Suppress("EXTENSION_SHADOWED_BY_MEMBER", "OVERRIDE_BY_INLINE", "NOTHING_TO_INLINE") @Suppress("EXTENSION_SHADOWED_BY_MEMBER", "OVERRIDE_BY_INLINE", "NOTHING_TO_INLINE")
public object FloatField : ExtendedField<Float>, Norm<Float, Float> { public object FloatField : ExtendedField<Float>, Norm<Float, Float> {
public override val zero: Float public override val zero: Float = 0.0f
get() = 0.0f public override val one: Float = 1.0f
public override val one: Float
get() = 1.0f
override fun number(value: Number): Float = value.toFloat() override fun number(value: Number): Float = value.toFloat()
@ -125,7 +120,7 @@ public object FloatField : ExtendedField<Float>, Norm<Float, Float> {
} }
public override inline fun add(a: Float, b: Float): Float = a + b public override inline fun add(a: Float, b: Float): Float = a + b
public override inline fun multiply(a: Float, k: Number): Float = a * k.toFloat() override fun scale(a: Float, value: Double): Float = a * value.toFloat()
public override inline fun multiply(a: Float, b: Float): Float = a * b public override inline fun multiply(a: Float, b: Float): Float = a * b
@ -170,12 +165,8 @@ public object IntRing : Ring<Int>, Norm<Int, Int>, NumericAlgebra<Int> {
get() = 1 get() = 1
override fun number(value: Number): Int = value.toInt() override fun number(value: Number): Int = value.toInt()
public override inline fun add(a: Int, b: Int): Int = a + b public override inline fun add(a: Int, b: Int): Int = a + b
public override inline fun multiply(a: Int, k: Number): Int = k.toInt() * a
public override inline fun multiply(a: Int, b: Int): Int = a * b public override inline fun multiply(a: Int, b: Int): Int = a * b
public override inline fun norm(arg: Int): Int = abs(arg) public override inline fun norm(arg: Int): Int = abs(arg)
public override inline fun Int.unaryMinus(): Int = -this public override inline fun Int.unaryMinus(): Int = -this
@ -196,12 +187,8 @@ public object ShortRing : Ring<Short>, Norm<Short, Short>, NumericAlgebra<Short>
get() = 1 get() = 1
override fun number(value: Number): Short = value.toShort() override fun number(value: Number): Short = value.toShort()
public override inline fun add(a: Short, b: Short): Short = (a + b).toShort() public override inline fun add(a: Short, b: Short): Short = (a + b).toShort()
public override inline fun multiply(a: Short, k: Number): Short = (a * k.toShort()).toShort()
public override inline fun multiply(a: Short, b: Short): Short = (a * b).toShort() public override inline fun multiply(a: Short, b: Short): Short = (a * b).toShort()
public override fun norm(arg: Short): Short = if (arg > 0) arg else (-arg).toShort() public override fun norm(arg: Short): Short = if (arg > 0) arg else (-arg).toShort()
public override inline fun Short.unaryMinus(): Short = (-this).toShort() public override inline fun Short.unaryMinus(): Short = (-this).toShort()
@ -222,12 +209,8 @@ public object ByteRing : Ring<Byte>, Norm<Byte, Byte>, NumericAlgebra<Byte> {
get() = 1 get() = 1
override fun number(value: Number): Byte = value.toByte() override fun number(value: Number): Byte = value.toByte()
public override inline fun add(a: Byte, b: Byte): Byte = (a + b).toByte() public override inline fun add(a: Byte, b: Byte): Byte = (a + b).toByte()
public override inline fun multiply(a: Byte, k: Number): Byte = (a * k.toByte()).toByte()
public override inline fun multiply(a: Byte, b: Byte): Byte = (a * b).toByte() public override inline fun multiply(a: Byte, b: Byte): Byte = (a * b).toByte()
public override fun norm(arg: Byte): Byte = if (arg > 0) arg else (-arg).toByte() public override fun norm(arg: Byte): Byte = if (arg > 0) arg else (-arg).toByte()
public override inline fun Byte.unaryMinus(): Byte = (-this).toByte() public override inline fun Byte.unaryMinus(): Byte = (-this).toByte()
@ -248,12 +231,8 @@ public object LongRing : Ring<Long>, Norm<Long, Long>, NumericAlgebra<Long> {
get() = 1L get() = 1L
override fun number(value: Number): Long = value.toLong() override fun number(value: Number): Long = value.toLong()
public override inline fun add(a: Long, b: Long): Long = a + b public override inline fun add(a: Long, b: Long): Long = a + b
public override inline fun multiply(a: Long, k: Number): Long = a * k.toLong()
public override inline fun multiply(a: Long, b: Long): Long = a * b public override inline fun multiply(a: Long, b: Long): Long = a * b
public override fun norm(arg: Long): Long = abs(arg) public override fun norm(arg: Long): Long = abs(arg)
public override inline fun Long.unaryMinus(): Long = (-this) public override inline fun Long.unaryMinus(): Long = (-this)

@ -100,9 +100,29 @@ public fun <T> Buffer<T>.asSequence(): Sequence<T> = Sequence(::iterator)
public fun <T> Buffer<T>.asIterable(): Iterable<T> = Iterable(::iterator) public fun <T> Buffer<T>.asIterable(): Iterable<T> = Iterable(::iterator)
/** /**
* Converts this [Buffer] to a new [List] * Returns a new [List] containing all elements of this buffer.
*/ */
public fun <T> Buffer<T>.toList(): List<T> = asSequence().toList() public fun <T> Buffer<T>.toList(): List<T> = when (this) {
is ArrayBuffer<T> -> array.toList()
is ListBuffer<T> -> list.toList()
is MutableListBuffer<T> -> list.toList()
else -> asSequence().toList()
}
/**
* Returns a new [MutableList] filled with all elements of this buffer.
*/
public fun <T> Buffer<T>.toMutableList(): MutableList<T> = when (this) {
is ArrayBuffer<T> -> array.toMutableList()
is ListBuffer<T> -> list.toMutableList()
is MutableListBuffer<T> -> list.toMutableList()
else -> asSequence().toMutableList()
}
/**
* Returns a new [Array] containing all elements of this buffer.
*/
public inline fun <reified T> Buffer<T>.toTypedArray(): Array<T> = asSequence().toList().toTypedArray()
/** /**
* Returns an [IntRange] of the valid indices for this [Buffer]. * Returns an [IntRange] of the valid indices for this [Buffer].
@ -222,7 +242,7 @@ public inline class MutableListBuffer<T>(public val list: MutableList<T>) : Muta
* @param T the type of elements contained in the buffer. * @param T the type of elements contained in the buffer.
* @property array The underlying array. * @property array The underlying array.
*/ */
public class ArrayBuffer<T>(private val array: Array<T>) : MutableBuffer<T> { public class ArrayBuffer<T>(internal val array: Array<T>) : MutableBuffer<T> {
// Can't inline because array is invariant // Can't inline because array is invariant
override val size: Int override val size: Int
get() = array.size get() = array.size

@ -13,10 +13,10 @@ internal class BufferAccessor2D<T : Any>(
public val colNum: Int, public val colNum: Int,
val factory: MutableBufferFactory<T>, val factory: MutableBufferFactory<T>,
) { ) {
public operator fun Buffer<T>.get(i: Int, j: Int): T = get(i + colNum * j) public operator fun Buffer<T>.get(i: Int, j: Int): T = get(i*colNum + j)
public operator fun MutableBuffer<T>.set(i: Int, j: Int, value: T) { public operator fun MutableBuffer<T>.set(i: Int, j: Int, value: T) {
set(i + colNum * j, value) set(i*colNum + j, value)
} }
public inline fun create(crossinline init: (i: Int, j: Int) -> T): MutableBuffer<T> = public inline fun create(crossinline init: (i: Int, j: Int) -> T): MutableBuffer<T> =
@ -25,7 +25,7 @@ internal class BufferAccessor2D<T : Any>(
public fun create(mat: Structure2D<T>): MutableBuffer<T> = create { i, j -> mat[i, j] } public fun create(mat: Structure2D<T>): MutableBuffer<T> = create { i, j -> mat[i, j] }
//TODO optimize wrapper //TODO optimize wrapper
public fun MutableBuffer<T>.collect(): Structure2D<T> = NDStructure.build( public fun MutableBuffer<T>.collect(): Structure2D<T> = NDStructure.buffered(
DefaultStrides(intArrayOf(rowNum, colNum)), DefaultStrides(intArrayOf(rowNum, colNum)),
factory factory
) { (i, j) -> ) { (i, j) ->

@ -36,10 +36,12 @@ public inline fun FloatBuffer(size: Int, init: (Int) -> Float): FloatBuffer = Fl
public fun FloatBuffer(vararg floats: Float): FloatBuffer = FloatBuffer(floats) public fun FloatBuffer(vararg floats: Float): FloatBuffer = FloatBuffer(floats)
/** /**
* Returns a [FloatArray] containing all of the elements of this [MutableBuffer]. * Returns a new [FloatArray] containing all of the elements of this [Buffer].
*/ */
public val MutableBuffer<out Float>.array: FloatArray public fun Buffer<Float>.toFloatArray(): FloatArray = when(this) {
get() = (if (this is FloatBuffer) array else FloatArray(size) { get(it) }) is FloatBuffer -> array.copyOf()
else -> FloatArray(size, ::get)
}
/** /**
* Returns [FloatBuffer] over this array. * Returns [FloatBuffer] over this array.

@ -35,10 +35,12 @@ public inline fun IntBuffer(size: Int, init: (Int) -> Int): IntBuffer = IntBuffe
public fun IntBuffer(vararg ints: Int): IntBuffer = IntBuffer(ints) public fun IntBuffer(vararg ints: Int): IntBuffer = IntBuffer(ints)
/** /**
* Returns a [IntArray] containing all of the elements of this [MutableBuffer]. * Returns a new [IntArray] containing all of the elements of this [Buffer].
*/ */
public val MutableBuffer<out Int>.array: IntArray public fun Buffer<Int>.toIntArray(): IntArray = when(this) {
get() = (if (this is IntBuffer) array else IntArray(size) { get(it) }) is IntBuffer -> array.copyOf()
else -> IntArray(size, ::get)
}
/** /**
* Returns [IntBuffer] over this array. * Returns [IntBuffer] over this array.

@ -35,10 +35,12 @@ public inline fun LongBuffer(size: Int, init: (Int) -> Long): LongBuffer = LongB
public fun LongBuffer(vararg longs: Long): LongBuffer = LongBuffer(longs) public fun LongBuffer(vararg longs: Long): LongBuffer = LongBuffer(longs)
/** /**
* Returns a [IntArray] containing all of the elements of this [MutableBuffer]. * Returns a new [LongArray] containing all of the elements of this [Buffer].
*/ */
public val MutableBuffer<out Long>.array: LongArray public fun Buffer<Long>.toLongArray(): LongArray = when(this) {
get() = (if (this is LongBuffer) array else LongArray(size) { get(it) }) is LongBuffer -> array.copyOf()
else -> LongArray(size, ::get)
}
/** /**
* Returns [LongBuffer] over this array. * Returns [LongBuffer] over this array.

@ -40,10 +40,12 @@ public fun RealBuffer(vararg doubles: Double): RealBuffer = RealBuffer(doubles)
public fun RealBuffer.contentEquals(vararg doubles: Double): Boolean = array.contentEquals(doubles) public fun RealBuffer.contentEquals(vararg doubles: Double): Boolean = array.contentEquals(doubles)
/** /**
* Returns a [DoubleArray] containing all of the elements of this [MutableBuffer]. * Returns a new [DoubleArray] containing all of the elements of this [Buffer].
*/ */
public val MutableBuffer<out Double>.array: DoubleArray public fun Buffer<Double>.toDoubleArray(): DoubleArray = when(this) {
get() = (if (this is RealBuffer) array else DoubleArray(size) { get(it) }) is RealBuffer -> array.copyOf()
else -> DoubleArray(size, ::get)
}
/** /**
* Returns [RealBuffer] over this array. * Returns [RealBuffer] over this array.

@ -8,6 +8,12 @@ import kotlin.math.*
* [ExtendedFieldOperations] over [RealBuffer]. * [ExtendedFieldOperations] over [RealBuffer].
*/ */
public object RealBufferFieldOperations : ExtendedFieldOperations<Buffer<Double>> { public object RealBufferFieldOperations : ExtendedFieldOperations<Buffer<Double>> {
override fun Buffer<Double>.unaryMinus(): RealBuffer = if (this is RealBuffer) {
RealBuffer(size) { -array[it] }
} else {
RealBuffer(size) { -get(it) }
}
public override fun add(a: Buffer<Double>, b: Buffer<Double>): RealBuffer { public override fun add(a: Buffer<Double>, b: Buffer<Double>): RealBuffer {
require(b.size == a.size) { require(b.size == a.size) {
"The size of the first buffer ${a.size} should be the same as for second one: ${b.size} " "The size of the first buffer ${a.size} should be the same as for second one: ${b.size} "
@ -19,15 +25,24 @@ public object RealBufferFieldOperations : ExtendedFieldOperations<Buffer<Double>
RealBuffer(DoubleArray(a.size) { aArray[it] + bArray[it] }) RealBuffer(DoubleArray(a.size) { aArray[it] + bArray[it] })
} else RealBuffer(DoubleArray(a.size) { a[it] + b[it] }) } else RealBuffer(DoubleArray(a.size) { a[it] + b[it] })
} }
//
public override fun multiply(a: Buffer<Double>, k: Number): RealBuffer { // public override fun multiply(a: Buffer<Double>, k: Number): RealBuffer {
val kValue = k.toDouble() // val kValue = k.toDouble()
//
return if (a is RealBuffer) { // return if (a is RealBuffer) {
val aArray = a.array // val aArray = a.array
RealBuffer(DoubleArray(a.size) { aArray[it] * kValue }) // RealBuffer(DoubleArray(a.size) { aArray[it] * kValue })
} else RealBuffer(DoubleArray(a.size) { a[it] * kValue }) // } else RealBuffer(DoubleArray(a.size) { a[it] * kValue })
} // }
//
// public override fun divide(a: Buffer<Double>, k: Number): RealBuffer {
// val kValue = k.toDouble()
//
// return if (a is RealBuffer) {
// val aArray = a.array
// RealBuffer(DoubleArray(a.size) { aArray[it] / kValue })
// } else RealBuffer(DoubleArray(a.size) { a[it] / kValue })
// }
public override fun multiply(a: Buffer<Double>, b: Buffer<Double>): RealBuffer { public override fun multiply(a: Buffer<Double>, b: Buffer<Double>): RealBuffer {
require(b.size == a.size) { require(b.size == a.size) {
@ -152,14 +167,22 @@ public class RealBufferField(public val size: Int) : ExtendedField<Buffer<Double
override fun number(value: Number): Buffer<Double> = RealBuffer(size) { value.toDouble() } override fun number(value: Number): Buffer<Double> = RealBuffer(size) { value.toDouble() }
override fun Buffer<Double>.unaryMinus(): Buffer<Double> = RealBufferFieldOperations.run {
-this@unaryMinus
}
public override fun add(a: Buffer<Double>, b: Buffer<Double>): RealBuffer { public override fun add(a: Buffer<Double>, b: Buffer<Double>): RealBuffer {
require(a.size == size) { "The buffer size ${a.size} does not match context size $size" } require(a.size == size) { "The buffer size ${a.size} does not match context size $size" }
return RealBufferFieldOperations.add(a, b) return RealBufferFieldOperations.add(a, b)
} }
public override fun multiply(a: Buffer<Double>, k: Number): RealBuffer { public override fun scale(a: Buffer<Double>, value: Double): RealBuffer {
require(a.size == size) { "The buffer size ${a.size} does not match context size $size" } require(a.size == size) { "The buffer size ${a.size} does not match context size $size" }
return RealBufferFieldOperations.multiply(a, k)
return if (a is RealBuffer) {
val aArray = a.array
RealBuffer(DoubleArray(a.size) { aArray[it] * value })
} else RealBuffer(DoubleArray(a.size) { a[it] * value })
} }
public override fun multiply(a: Buffer<Double>, b: Buffer<Double>): RealBuffer { public override fun multiply(a: Buffer<Double>, b: Buffer<Double>): RealBuffer {

@ -33,10 +33,12 @@ public inline fun ShortBuffer(size: Int, init: (Int) -> Short): ShortBuffer = Sh
public fun ShortBuffer(vararg shorts: Short): ShortBuffer = ShortBuffer(shorts) public fun ShortBuffer(vararg shorts: Short): ShortBuffer = ShortBuffer(shorts)
/** /**
* Returns a [ShortArray] containing all of the elements of this [MutableBuffer]. * Returns a new [ShortArray] containing all of the elements of this [Buffer].
*/ */
public val MutableBuffer<out Short>.array: ShortArray public fun Buffer<Short>.toShortArray(): ShortArray = when(this) {
get() = (if (this is ShortBuffer) array else ShortArray(size) { get(it) }) is ShortBuffer -> array.copyOf()
else -> ShortArray(size, ::get)
}
/** /**
* Returns [ShortBuffer] over this array. * Returns [ShortBuffer] over this array.

@ -0,0 +1,7 @@
package space.kscience.kmath.structures
public inline fun <T : Any, reified R : Any> Buffer<T>.map(
bufferFactory: BufferFactory<R> = Buffer.Companion::auto,
crossinline block: (T) -> R,
): Buffer<R> = bufferFactory(size) { block(get(it)) }

@ -11,9 +11,7 @@ class ExpressionFieldTest {
@Test @Test
fun testExpression() { fun testExpression() {
val context = FunctionalExpressionField(RealField) val expression = FunctionalExpressionField(RealField).invoke {
val expression = context {
val x by binding() val x by binding()
x * x + 2 * x + one x * x + 2 * x + one
} }

@ -1,23 +1,24 @@
package space.kscience.kmath.linear package space.kscience.kmath.linear
import space.kscience.kmath.misc.UnstableKMathAPI
import space.kscience.kmath.nd.NDStructure import space.kscience.kmath.nd.NDStructure
import space.kscience.kmath.nd.as2D import space.kscience.kmath.nd.as2D
import space.kscience.kmath.operations.invoke
import kotlin.test.Test import kotlin.test.Test
import kotlin.test.assertEquals import kotlin.test.assertEquals
@UnstableKMathAPI
@Suppress("UNUSED_VARIABLE") @Suppress("UNUSED_VARIABLE")
class MatrixTest { class MatrixTest {
@Test @Test
fun testTranspose() { fun testTranspose() {
val matrix = MatrixContext.real.one(3, 3) val matrix = LinearSpace.real.one(3, 3)
val transposed = matrix.transpose() val transposed = matrix.transpose()
assertEquals(matrix, transposed) assertEquals(matrix, transposed)
} }
@Test @Test
fun testBuilder() { fun testBuilder() {
val matrix = Matrix.build(2, 3)( val matrix = LinearSpace.real.matrix(2, 3)(
1.0, 0.0, 0.0, 1.0, 0.0, 0.0,
0.0, 1.0, 2.0 0.0, 1.0, 2.0
) )
@ -39,7 +40,7 @@ class MatrixTest {
infix fun Matrix<Double>.pow(power: Int): Matrix<Double> { infix fun Matrix<Double>.pow(power: Int): Matrix<Double> {
var res = this var res = this
repeat(power - 1) { repeat(power - 1) {
res = RealMatrixContext.invoke { res dot this@pow } res = LinearSpace.real.run { res dot this@pow }
} }
return res return res
} }
@ -52,7 +53,7 @@ class MatrixTest {
val firstMatrix = NDStructure.auto(2, 3) { (i, j) -> (i + j).toDouble() }.as2D() val firstMatrix = NDStructure.auto(2, 3) { (i, j) -> (i + j).toDouble() }.as2D()
val secondMatrix = NDStructure.auto(3, 2) { (i, j) -> (i + j).toDouble() }.as2D() val secondMatrix = NDStructure.auto(3, 2) { (i, j) -> (i + j).toDouble() }.as2D()
MatrixContext.real.run { LinearSpace.real.run {
// val firstMatrix = produce(2, 3) { i, j -> (i + j).toDouble() } // val firstMatrix = produce(2, 3) { i, j -> (i + j).toDouble() }
// val secondMatrix = produce(3, 2) { i, j -> (i + j).toDouble() } // val secondMatrix = produce(3, 2) { i, j -> (i + j).toDouble() }
val result = firstMatrix dot secondMatrix val result = firstMatrix dot secondMatrix

@ -1,29 +1,31 @@
package space.kscience.kmath.linear package space.kscience.kmath.linear
import space.kscience.kmath.misc.UnstableKMathAPI
import kotlin.test.Test import kotlin.test.Test
import kotlin.test.assertEquals import kotlin.test.assertEquals
@UnstableKMathAPI
class RealLUSolverTest { class RealLUSolverTest {
@Test @Test
fun testInvertOne() { fun testInvertOne() {
val matrix = MatrixContext.real.one(2, 2) val matrix = LinearSpace.real.one(2, 2)
val inverted = MatrixContext.real.inverseWithLup(matrix) val inverted = LinearSpace.real.inverseWithLup(matrix)
assertEquals(matrix, inverted) assertEquals(matrix, inverted)
} }
@Test @Test
fun testDecomposition() { fun testDecomposition() {
val matrix = Matrix.square( LinearSpace.real.run {
3.0, 1.0, val matrix = matrix(2,2)(
1.0, 3.0 3.0, 1.0,
) 2.0, 3.0
)
MatrixContext.real.run {
val lup = lup(matrix) val lup = lup(matrix)
//Check determinant //Check determinant
assertEquals(8.0, lup.determinant) assertEquals(7.0, lup.determinant)
assertEquals(lup.p dot matrix, lup.l dot lup.u) assertEquals(lup.p dot matrix, lup.l dot lup.u)
} }
@ -31,14 +33,14 @@ class RealLUSolverTest {
@Test @Test
fun testInvert() { fun testInvert() {
val matrix = Matrix.square( val matrix = LinearSpace.real.matrix(2,2)(
3.0, 1.0, 3.0, 1.0,
1.0, 3.0 1.0, 3.0
) )
val inverted = MatrixContext.real.inverseWithLup(matrix) val inverted = LinearSpace.real.inverseWithLup(matrix)
val expected = Matrix.square( val expected = LinearSpace.real.matrix(2,2)(
0.375, -0.125, 0.375, -0.125,
-0.125, 0.375 -0.125, 0.375
) )

@ -1,5 +1,6 @@
package space.kscience.kmath.structures package space.kscience.kmath.structures
import space.kscience.kmath.linear.LinearSpace
import space.kscience.kmath.nd.* import space.kscience.kmath.nd.*
import space.kscience.kmath.operations.Norm import space.kscience.kmath.operations.Norm
import space.kscience.kmath.operations.invoke import space.kscience.kmath.operations.invoke
@ -33,7 +34,9 @@ class NumberNDFieldTest {
@Test @Test
fun testGeneration() { fun testGeneration() {
val array = Structure2D.real(3, 3) { i, j -> (i * 10 + j).toDouble() } val array = LinearSpace.real.buildMatrix(3, 3) { i, j ->
(i * 10 + j).toDouble()
}
for (i in 0..2) { for (i in 0..2) {
for (j in 0..2) { for (j in 0..2) {

@ -5,8 +5,9 @@ import space.kscience.kmath.operations.invoke
import kotlin.test.assertEquals import kotlin.test.assertEquals
import kotlin.test.assertNotEquals import kotlin.test.assertNotEquals
internal class FieldVerifier<T>(override val algebra: Field<T>, a: T, b: T, c: T, x: Number) : internal class FieldVerifier<T, A : Field<T>>(
RingVerifier<T>(algebra, a, b, c, x) { algebra: A, a: T, b: T, c: T, x: Number,
) : RingVerifier<T, A>(algebra, a, b, c, x) {
override fun verify() { override fun verify() {
super.verify() super.verify()

@ -1,11 +1,13 @@
package space.kscience.kmath.testutils package space.kscience.kmath.testutils
import space.kscience.kmath.operations.Ring import space.kscience.kmath.operations.Ring
import space.kscience.kmath.operations.ScaleOperations
import space.kscience.kmath.operations.invoke import space.kscience.kmath.operations.invoke
import kotlin.test.assertEquals import kotlin.test.assertEquals
internal open class RingVerifier<T>(override val algebra: Ring<T>, a: T, b: T, c: T, x: Number) : internal open class RingVerifier<T, A>(algebra: A, a: T, b: T, c: T, x: Number) :
SpaceVerifier<T>(algebra, a, b, c, x) { SpaceVerifier<T, A>(algebra, a, b, c, x) where A : Ring<T>, A : ScaleOperations<T> {
override fun verify() { override fun verify() {
super.verify() super.verify()

@ -1,18 +1,18 @@
package space.kscience.kmath.testutils package space.kscience.kmath.testutils
import space.kscience.kmath.operations.Space import space.kscience.kmath.operations.Group
import space.kscience.kmath.operations.ScaleOperations
import space.kscience.kmath.operations.invoke import space.kscience.kmath.operations.invoke
import kotlin.test.assertEquals import kotlin.test.assertEquals
import kotlin.test.assertNotEquals import kotlin.test.assertNotEquals
internal open class SpaceVerifier<T>( internal open class SpaceVerifier<T, out S>(
override val algebra: Space<T>, override val algebra: S,
val a: T, val a: T,
val b: T, val b: T,
val c: T, val c: T,
val x: Number val x: Number,
) : ) : AlgebraicVerifier<T, Group<T>> where S : Group<T>, S : ScaleOperations<T> {
AlgebraicVerifier<T, Space<T>> {
override fun verify() { override fun verify() {
algebra { algebra {
assertEquals(a + b + c, a + (b + c), "Addition in $algebra is not associative.") assertEquals(a + b + c, a + (b + c), "Addition in $algebra is not associative.")

@ -7,19 +7,16 @@ import java.math.MathContext
/** /**
* A field over [BigInteger]. * A field over [BigInteger].
*/ */
public object JBigIntegerField : Field<BigInteger>, NumericAlgebra<BigInteger> { public object JBigIntegerField : Ring<BigInteger>, NumericAlgebra<BigInteger> {
public override val zero: BigInteger public override val zero: BigInteger get() = BigInteger.ZERO
get() = BigInteger.ZERO
public override val one: BigInteger public override val one: BigInteger get() = BigInteger.ONE
get() = BigInteger.ONE
public override fun number(value: Number): BigInteger = BigInteger.valueOf(value.toLong()) public override fun number(value: Number): BigInteger = BigInteger.valueOf(value.toLong())
public override fun divide(a: BigInteger, b: BigInteger): BigInteger = a.div(b)
public override fun add(a: BigInteger, b: BigInteger): BigInteger = a.add(b) public override fun add(a: BigInteger, b: BigInteger): BigInteger = a.add(b)
public override operator fun BigInteger.minus(b: BigInteger): BigInteger = subtract(b) public override operator fun BigInteger.minus(b: BigInteger): BigInteger = subtract(b)
public override fun multiply(a: BigInteger, k: Number): BigInteger = a.multiply(k.toInt().toBigInteger())
public override fun multiply(a: BigInteger, b: BigInteger): BigInteger = a.multiply(b) public override fun multiply(a: BigInteger, b: BigInteger): BigInteger = a.multiply(b)
public override operator fun BigInteger.unaryMinus(): BigInteger = negate() public override operator fun BigInteger.unaryMinus(): BigInteger = negate()
} }
@ -30,7 +27,7 @@ public object JBigIntegerField : Field<BigInteger>, NumericAlgebra<BigInteger> {
*/ */
public abstract class JBigDecimalFieldBase internal constructor( public abstract class JBigDecimalFieldBase internal constructor(
private val mathContext: MathContext = MathContext.DECIMAL64, private val mathContext: MathContext = MathContext.DECIMAL64,
) : Field<BigDecimal>, PowerOperations<BigDecimal>, NumericAlgebra<BigDecimal> { ) : Field<BigDecimal>, PowerOperations<BigDecimal>, NumericAlgebra<BigDecimal>, ScaleOperations<BigDecimal> {
public override val zero: BigDecimal public override val zero: BigDecimal
get() = BigDecimal.ZERO get() = BigDecimal.ZERO
@ -41,8 +38,8 @@ public abstract class JBigDecimalFieldBase internal constructor(
public override operator fun BigDecimal.minus(b: BigDecimal): BigDecimal = subtract(b) public override operator fun BigDecimal.minus(b: BigDecimal): BigDecimal = subtract(b)
public override fun number(value: Number): BigDecimal = BigDecimal.valueOf(value.toDouble()) public override fun number(value: Number): BigDecimal = BigDecimal.valueOf(value.toDouble())
public override fun multiply(a: BigDecimal, k: Number): BigDecimal = public override fun scale(a: BigDecimal, value: Double): BigDecimal =
a.multiply(k.toDouble().toBigDecimal(mathContext), mathContext) a.multiply(value.toBigDecimal(mathContext), mathContext)
public override fun multiply(a: BigDecimal, b: BigDecimal): BigDecimal = a.multiply(b, mathContext) public override fun multiply(a: BigDecimal, b: BigDecimal): BigDecimal = a.multiply(b, mathContext)
public override fun divide(a: BigDecimal, b: BigDecimal): BigDecimal = a.divide(b, mathContext) public override fun divide(a: BigDecimal, b: BigDecimal): BigDecimal = a.divide(b, mathContext)

@ -5,16 +5,16 @@ import kotlinx.coroutines.flow.Flow
import kotlinx.coroutines.flow.map import kotlinx.coroutines.flow.map
import kotlinx.coroutines.flow.runningReduce import kotlinx.coroutines.flow.runningReduce
import kotlinx.coroutines.flow.scan import kotlinx.coroutines.flow.scan
import space.kscience.kmath.operations.Space import space.kscience.kmath.operations.Group
import space.kscience.kmath.operations.SpaceOperations import space.kscience.kmath.operations.GroupOperations
import space.kscience.kmath.operations.ScaleOperations
import space.kscience.kmath.operations.invoke import space.kscience.kmath.operations.invoke
@ExperimentalCoroutinesApi public fun <T> Flow<T>.cumulativeSum(group: GroupOperations<T>): Flow<T> =
public fun <T> Flow<T>.cumulativeSum(space: SpaceOperations<T>): Flow<T> = group { runningReduce { sum, element -> sum + element } }
space { runningReduce { sum, element -> sum + element } }
@ExperimentalCoroutinesApi @ExperimentalCoroutinesApi
public fun <T> Flow<T>.mean(space: Space<T>): Flow<T> = space { public fun <T, S> Flow<T>.mean(algebra: S): Flow<T> where S : Group<T>, S : ScaleOperations<T> = algebra {
data class Accumulator(var sum: T, var num: Int) data class Accumulator(var sum: T, var num: Int)
scan(Accumulator(zero, 0)) { sum, element -> scan(Accumulator(zero, 0)) { sum, element ->

@ -2,7 +2,8 @@ package space.kscience.kmath.dimensions
import space.kscience.kmath.linear.* import space.kscience.kmath.linear.*
import space.kscience.kmath.nd.Structure2D import space.kscience.kmath.nd.Structure2D
import space.kscience.kmath.operations.invoke import space.kscience.kmath.operations.RealField
import space.kscience.kmath.operations.Ring
/** /**
* A matrix with compile-time controlled dimension * A matrix with compile-time controlled dimension
@ -77,7 +78,7 @@ public inline class DPointWrapper<T, D : Dimension>(public val point: Point<T>)
/** /**
* Basic operations on dimension-safe matrices. Operates on [Matrix] * Basic operations on dimension-safe matrices. Operates on [Matrix]
*/ */
public inline class DMatrixContext<T : Any>(public val context: MatrixContext<T, Matrix<T>>) { public inline class DMatrixContext<T : Any, out A : Ring<T>>(public val context: LinearSpace<T, A>) {
public inline fun <reified R : Dimension, reified C : Dimension> Matrix<T>.coerce(): DMatrix<T, R, C> { public inline fun <reified R : Dimension, reified C : Dimension> Matrix<T>.coerce(): DMatrix<T, R, C> {
require(rowNum == Dimension.dim<R>().toInt()) { require(rowNum == Dimension.dim<R>().toInt()) {
"Row number mismatch: expected ${Dimension.dim<R>()} but found $rowNum" "Row number mismatch: expected ${Dimension.dim<R>()} but found $rowNum"
@ -93,17 +94,19 @@ public inline class DMatrixContext<T : Any>(public val context: MatrixContext<T,
/** /**
* Produce a matrix with this context and given dimensions * Produce a matrix with this context and given dimensions
*/ */
public inline fun <reified R : Dimension, reified C : Dimension> produce(noinline initializer: (i: Int, j: Int) -> T): DMatrix<T, R, C> { public inline fun <reified R : Dimension, reified C : Dimension> produce(
noinline initializer: A.(i: Int, j: Int) -> T
): DMatrix<T, R, C> {
val rows = Dimension.dim<R>() val rows = Dimension.dim<R>()
val cols = Dimension.dim<C>() val cols = Dimension.dim<C>()
return context.produce(rows.toInt(), cols.toInt(), initializer).coerce<R, C>() return context.buildMatrix(rows.toInt(), cols.toInt(), initializer).coerce<R, C>()
} }
public inline fun <reified D : Dimension> point(noinline initializer: (Int) -> T): DPoint<T, D> { public inline fun <reified D : Dimension> point(noinline initializer: A.(Int) -> T): DPoint<T, D> {
val size = Dimension.dim<D>() val size = Dimension.dim<D>()
return DPoint.coerceUnsafe( return DPoint.coerceUnsafe(
context.point( context.buildVector(
size.toInt(), size.toInt(),
initializer initializer
) )
@ -112,31 +115,31 @@ public inline class DMatrixContext<T : Any>(public val context: MatrixContext<T,
public inline infix fun <reified R1 : Dimension, reified C1 : Dimension, reified C2 : Dimension> DMatrix<T, R1, C1>.dot( public inline infix fun <reified R1 : Dimension, reified C1 : Dimension, reified C2 : Dimension> DMatrix<T, R1, C1>.dot(
other: DMatrix<T, C1, C2>, other: DMatrix<T, C1, C2>,
): DMatrix<T, R1, C2> = context { this@dot dot other }.coerce() ): DMatrix<T, R1, C2> = context.run { this@dot dot other }.coerce()
public inline infix fun <reified R : Dimension, reified C : Dimension> DMatrix<T, R, C>.dot(vector: DPoint<T, C>): DPoint<T, R> = public inline infix fun <reified R : Dimension, reified C : Dimension> DMatrix<T, R, C>.dot(vector: DPoint<T, C>): DPoint<T, R> =
DPoint.coerceUnsafe(context { this@dot dot vector }) DPoint.coerceUnsafe(context.run { this@dot dot vector })
public inline operator fun <reified R : Dimension, reified C : Dimension> DMatrix<T, R, C>.times(value: T): DMatrix<T, R, C> = public inline operator fun <reified R : Dimension, reified C : Dimension> DMatrix<T, R, C>.times(value: T): DMatrix<T, R, C> =
context { this@times.times(value) }.coerce() context.run { this@times.times(value) }.coerce()
public inline operator fun <reified R : Dimension, reified C : Dimension> T.times(m: DMatrix<T, R, C>): DMatrix<T, R, C> = public inline operator fun <reified R : Dimension, reified C : Dimension> T.times(m: DMatrix<T, R, C>): DMatrix<T, R, C> =
m * this m * this
public inline operator fun <reified R : Dimension, reified C : Dimension> DMatrix<T, C, R>.plus(other: DMatrix<T, C, R>): DMatrix<T, C, R> = public inline operator fun <reified R : Dimension, reified C : Dimension> DMatrix<T, C, R>.plus(other: DMatrix<T, C, R>): DMatrix<T, C, R> =
context { this@plus + other }.coerce() context.run { this@plus + other }.coerce()
public inline operator fun <reified R : Dimension, reified C : Dimension> DMatrix<T, C, R>.minus(other: DMatrix<T, C, R>): DMatrix<T, C, R> = public inline operator fun <reified R : Dimension, reified C : Dimension> DMatrix<T, C, R>.minus(other: DMatrix<T, C, R>): DMatrix<T, C, R> =
context { this@minus + other }.coerce() context.run { this@minus + other }.coerce()
public inline operator fun <reified R : Dimension, reified C : Dimension> DMatrix<T, C, R>.unaryMinus(): DMatrix<T, C, R> = public inline operator fun <reified R : Dimension, reified C : Dimension> DMatrix<T, C, R>.unaryMinus(): DMatrix<T, C, R> =
context { this@unaryMinus.unaryMinus() }.coerce() context.run { this@unaryMinus.unaryMinus() }.coerce()
public inline fun <reified R : Dimension, reified C : Dimension> DMatrix<T, C, R>.transpose(): DMatrix<T, R, C> = public inline fun <reified R : Dimension, reified C : Dimension> DMatrix<T, C, R>.transpose(): DMatrix<T, R, C> =
context { (this@transpose as Matrix<T>).transpose() }.coerce() context.run { (this@transpose as Matrix<T>).transpose() }.coerce()
public companion object { public companion object {
public val real: DMatrixContext<Double> = DMatrixContext(MatrixContext.real) public val real: DMatrixContext<Double, RealField> = DMatrixContext(LinearSpace.real)
} }
} }
@ -144,11 +147,11 @@ public inline class DMatrixContext<T : Any>(public val context: MatrixContext<T,
/** /**
* A square unit matrix * A square unit matrix
*/ */
public inline fun <reified D : Dimension> DMatrixContext<Double>.one(): DMatrix<Double, D, D> = produce { i, j -> public inline fun <reified D : Dimension> DMatrixContext<Double, RealField>.one(): DMatrix<Double, D, D> = produce { i, j ->
if (i == j) 1.0 else 0.0 if (i == j) 1.0 else 0.0
} }
public inline fun <reified R : Dimension, reified C : Dimension> DMatrixContext<Double>.zero(): DMatrix<Double, R, C> = public inline fun <reified R : Dimension, reified C : Dimension> DMatrixContext<Double, RealField>.zero(): DMatrix<Double, R, C> =
produce { _, _ -> produce { _, _ ->
0.0 0.0
} }

@ -0,0 +1,113 @@
package space.kscience.kmath.ejml
import org.ejml.simple.SimpleMatrix
import space.kscience.kmath.linear.*
import space.kscience.kmath.misc.UnstableKMathAPI
import space.kscience.kmath.nd.getFeature
import space.kscience.kmath.operations.RealField
/**
* Represents context of basic operations operating with [EjmlMatrix].
*
* @author Iaroslav Postovalov
*/
public object EjmlLinearSpace : LinearSpace<Double, RealField> {
override val elementAlgebra: RealField get() = RealField
/**
* Converts this matrix to EJML one.
*/
@OptIn(UnstableKMathAPI::class)
public fun Matrix<Double>.toEjml(): EjmlMatrix = when (val matrix = origin) {
is EjmlMatrix -> matrix
else -> buildMatrix(rowNum, colNum) { i, j -> get(i, j) }
}
/**
* Converts this vector to EJML one.
*/
public fun Point<Double>.toEjml(): EjmlVector = when (this) {
is EjmlVector -> this
else -> EjmlVector(SimpleMatrix(size, 1).also {
(0 until it.numRows()).forEach { row -> it[row, 0] = get(row) }
})
}
override fun buildMatrix(rows: Int, columns: Int, initializer: RealField.(i: Int, j: Int) -> Double): EjmlMatrix =
EjmlMatrix(SimpleMatrix(rows, columns).also {
(0 until rows).forEach { row ->
(0 until columns).forEach { col -> it[row, col] = RealField.initializer(row, col) }
}
})
override fun buildVector(size: Int, initializer: RealField.(Int) -> Double): Point<Double> =
EjmlVector(SimpleMatrix(size, 1).also {
(0 until it.numRows()).forEach { row -> it[row, 0] = RealField.initializer(row) }
})
private fun SimpleMatrix.wrapMatrix() = EjmlMatrix(this)
private fun SimpleMatrix.wrapVector() = EjmlVector(this)
override fun Matrix<Double>.unaryMinus(): Matrix<Double> = this * (-1.0)
public override fun Matrix<Double>.dot(other: Matrix<Double>): EjmlMatrix =
EjmlMatrix(toEjml().origin.mult(other.toEjml().origin))
public override fun Matrix<Double>.dot(vector: Point<Double>): EjmlVector =
EjmlVector(toEjml().origin.mult(vector.toEjml().origin))
public override operator fun Matrix<Double>.minus(other: Matrix<Double>): EjmlMatrix =
(toEjml().origin - other.toEjml().origin).wrapMatrix()
public override operator fun Matrix<Double>.times(value: Double): EjmlMatrix =
toEjml().origin.scale(value).wrapMatrix()
override fun Point<Double>.unaryMinus(): EjmlVector =
toEjml().origin.negative().wrapVector()
override fun Matrix<Double>.plus(other: Matrix<Double>): EjmlMatrix =
(toEjml().origin + other.toEjml().origin).wrapMatrix()
override fun Point<Double>.plus(other: Point<Double>): EjmlVector =
(toEjml().origin + other.toEjml().origin).wrapVector()
override fun Point<Double>.minus(other: Point<Double>): EjmlVector =
(toEjml().origin - other.toEjml().origin).wrapVector()
override fun Double.times(m: Matrix<Double>): EjmlMatrix =
m.toEjml().origin.scale(this).wrapMatrix()
override fun Point<Double>.times(value: Double): EjmlVector =
toEjml().origin.scale(value).wrapVector()
override fun Double.times(v: Point<Double>): EjmlVector =
v.toEjml().origin.scale(this).wrapVector()
}
/**
* Solves for X in the following equation: x = a^-1*b, where 'a' is base matrix and 'b' is an n by p matrix.
*
* @param a the base matrix.
* @param b n by p matrix.
* @return the solution for 'x' that is n by p.
* @author Iaroslav Postovalov
*/
public fun EjmlLinearSpace.solve(a: Matrix<Double>, b: Matrix<Double>): EjmlMatrix =
EjmlMatrix(a.toEjml().origin.solve(b.toEjml().origin))
/**
* Solves for X in the following equation: x = a^(-1)*b, where 'a' is base matrix and 'b' is an n by p matrix.
*
* @param a the base matrix.
* @param b n by p vector.
* @return the solution for 'x' that is n by p.
* @author Iaroslav Postovalov
*/
public fun EjmlLinearSpace.solve(a: Matrix<Double>, b: Point<Double>): EjmlVector =
EjmlVector(a.toEjml().origin.solve(b.toEjml().origin))
@OptIn(UnstableKMathAPI::class)
public fun EjmlMatrix.inverted(): EjmlMatrix = getFeature<InverseMatrixFeature<Double>>()!!.inverse as EjmlMatrix
public fun EjmlLinearSpace.inverse(matrix: Matrix<Double>): Matrix<Double> = matrix.toEjml().inverted()

@ -85,7 +85,7 @@ public class EjmlMatrix(public val origin: SimpleMatrix) : Matrix<Double> {
override fun equals(other: Any?): Boolean { override fun equals(other: Any?): Boolean {
if (this === other) return true if (this === other) return true
if (other !is Matrix<*>) return false if (other !is NDStructure<*>) return false
return NDStructure.contentEquals(this, other) return NDStructure.contentEquals(this, other)
} }

@ -1,88 +0,0 @@
package space.kscience.kmath.ejml
import org.ejml.simple.SimpleMatrix
import space.kscience.kmath.linear.*
import space.kscience.kmath.misc.UnstableKMathAPI
import space.kscience.kmath.nd.getFeature
/**
* Represents context of basic operations operating with [EjmlMatrix].
*
* @author Iaroslav Postovalov
*/
public object EjmlMatrixContext : MatrixContext<Double, EjmlMatrix> {
/**
* Converts this matrix to EJML one.
*/
@OptIn(UnstableKMathAPI::class)
public fun Matrix<Double>.toEjml(): EjmlMatrix = when (val matrix = origin) {
is EjmlMatrix -> matrix
else -> produce(rowNum, colNum) { i, j -> get(i, j) }
}
/**
* Converts this vector to EJML one.
*/
public fun Point<Double>.toEjml(): EjmlVector =
if (this is EjmlVector) this else EjmlVector(SimpleMatrix(size, 1).also {
(0 until it.numRows()).forEach { row -> it[row, 0] = get(row) }
})
override fun produce(rows: Int, columns: Int, initializer: (i: Int, j: Int) -> Double): EjmlMatrix =
EjmlMatrix(SimpleMatrix(rows, columns).also {
(0 until rows).forEach { row ->
(0 until columns).forEach { col -> it[row, col] = initializer(row, col) }
}
})
override fun point(size: Int, initializer: (Int) -> Double): Point<Double> =
EjmlVector(SimpleMatrix(size, 1).also {
(0 until it.numRows()).forEach { row -> it[row, 0] = initializer(row) }
})
public override fun Matrix<Double>.dot(other: Matrix<Double>): EjmlMatrix =
EjmlMatrix(toEjml().origin.mult(other.toEjml().origin))
public override fun Matrix<Double>.dot(vector: Point<Double>): EjmlVector =
EjmlVector(toEjml().origin.mult(vector.toEjml().origin))
public override fun add(a: Matrix<Double>, b: Matrix<Double>): EjmlMatrix =
EjmlMatrix(a.toEjml().origin + b.toEjml().origin)
public override operator fun Matrix<Double>.minus(b: Matrix<Double>): EjmlMatrix =
EjmlMatrix(toEjml().origin - b.toEjml().origin)
public override fun multiply(a: Matrix<Double>, k: Number): EjmlMatrix =
produce(a.rowNum, a.colNum) { i, j -> a[i, j] * k.toDouble() }
public override operator fun Matrix<Double>.times(value: Double): EjmlMatrix =
EjmlMatrix(toEjml().origin.scale(value))
}
/**
* Solves for X in the following equation: x = a^-1*b, where 'a' is base matrix and 'b' is an n by p matrix.
*
* @param a the base matrix.
* @param b n by p matrix.
* @return the solution for 'x' that is n by p.
* @author Iaroslav Postovalov
*/
public fun EjmlMatrixContext.solve(a: Matrix<Double>, b: Matrix<Double>): EjmlMatrix =
EjmlMatrix(a.toEjml().origin.solve(b.toEjml().origin))
/**
* Solves for X in the following equation: x = a^(-1)*b, where 'a' is base matrix and 'b' is an n by p matrix.
*
* @param a the base matrix.
* @param b n by p vector.
* @return the solution for 'x' that is n by p.
* @author Iaroslav Postovalov
*/
public fun EjmlMatrixContext.solve(a: Matrix<Double>, b: Point<Double>): EjmlVector =
EjmlVector(a.toEjml().origin.solve(b.toEjml().origin))
@OptIn(UnstableKMathAPI::class)
public fun EjmlMatrix.inverted(): EjmlMatrix = getFeature<InverseMatrixFeature<Double>>()!!.inverse as EjmlMatrix
public fun EjmlMatrixContext.inverse(matrix: Matrix<Double>): Matrix<Double> = matrix.toEjml().inverted()

@ -2,6 +2,7 @@ package space.kscience.kmath.real
import space.kscience.kmath.linear.* import space.kscience.kmath.linear.*
import space.kscience.kmath.misc.UnstableKMathAPI import space.kscience.kmath.misc.UnstableKMathAPI
import space.kscience.kmath.operations.RealField
import space.kscience.kmath.structures.Buffer import space.kscience.kmath.structures.Buffer
import space.kscience.kmath.structures.RealBuffer import space.kscience.kmath.structures.RealBuffer
import space.kscience.kmath.structures.asIterable import space.kscience.kmath.structures.asIterable
@ -21,15 +22,19 @@ import kotlin.math.pow
public typealias RealMatrix = Matrix<Double> public typealias RealMatrix = Matrix<Double>
public fun realMatrix(rowNum: Int, colNum: Int, initializer: (i: Int, j: Int) -> Double): RealMatrix = public fun realMatrix(rowNum: Int, colNum: Int, initializer: RealField.(i: Int, j: Int) -> Double): RealMatrix =
MatrixContext.real.produce(rowNum, colNum, initializer) LinearSpace.real.buildMatrix(rowNum, colNum, initializer)
@OptIn(UnstableKMathAPI::class)
public fun realMatrix(rowNum: Int, colNum: Int): MatrixBuilder<Double, RealField> =
LinearSpace.real.matrix(rowNum, colNum)
public fun Array<DoubleArray>.toMatrix(): RealMatrix { public fun Array<DoubleArray>.toMatrix(): RealMatrix {
return MatrixContext.real.produce(size, this[0].size) { row, col -> this[row][col] } return LinearSpace.real.buildMatrix(size, this[0].size) { row, col -> this@toMatrix[row][col] }
} }
public fun Sequence<DoubleArray>.toMatrix(): RealMatrix = toList().let { public fun Sequence<DoubleArray>.toMatrix(): RealMatrix = toList().let {
MatrixContext.real.produce(it.size, it[0].size) { row, col -> it[row][col] } LinearSpace.real.buildMatrix(it.size, it[0].size) { row, col -> it[row][col] }
} }
public fun RealMatrix.repeatStackVertical(n: Int): RealMatrix = public fun RealMatrix.repeatStackVertical(n: Int): RealMatrix =
@ -42,38 +47,38 @@ public fun RealMatrix.repeatStackVertical(n: Int): RealMatrix =
*/ */
public operator fun RealMatrix.times(double: Double): RealMatrix = public operator fun RealMatrix.times(double: Double): RealMatrix =
MatrixContext.real.produce(rowNum, colNum) { row, col -> LinearSpace.real.buildMatrix(rowNum, colNum) { row, col ->
this[row, col] * double get(row, col) * double
} }
public operator fun RealMatrix.plus(double: Double): RealMatrix = public operator fun RealMatrix.plus(double: Double): RealMatrix =
MatrixContext.real.produce(rowNum, colNum) { row, col -> LinearSpace.real.buildMatrix(rowNum, colNum) { row, col ->
this[row, col] + double get(row, col) + double
} }
public operator fun RealMatrix.minus(double: Double): RealMatrix = public operator fun RealMatrix.minus(double: Double): RealMatrix =
MatrixContext.real.produce(rowNum, colNum) { row, col -> LinearSpace.real.buildMatrix(rowNum, colNum) { row, col ->
this[row, col] - double get(row, col) - double
} }
public operator fun RealMatrix.div(double: Double): RealMatrix = public operator fun RealMatrix.div(double: Double): RealMatrix =
MatrixContext.real.produce(rowNum, colNum) { row, col -> LinearSpace.real.buildMatrix(rowNum, colNum) { row, col ->
this[row, col] / double get(row, col) / double
} }
public operator fun Double.times(matrix: RealMatrix): RealMatrix = public operator fun Double.times(matrix: RealMatrix): RealMatrix =
MatrixContext.real.produce(matrix.rowNum, matrix.colNum) { row, col -> LinearSpace.real.buildMatrix(matrix.rowNum, matrix.colNum) { row, col ->
this * matrix[row, col] this@times * matrix[row, col]
} }
public operator fun Double.plus(matrix: RealMatrix): RealMatrix = public operator fun Double.plus(matrix: RealMatrix): RealMatrix =
MatrixContext.real.produce(matrix.rowNum, matrix.colNum) { row, col -> LinearSpace.real.buildMatrix(matrix.rowNum, matrix.colNum) { row, col ->
this + matrix[row, col] this@plus + matrix[row, col]
} }
public operator fun Double.minus(matrix: RealMatrix): RealMatrix = public operator fun Double.minus(matrix: RealMatrix): RealMatrix =
MatrixContext.real.produce(matrix.rowNum, matrix.colNum) { row, col -> LinearSpace.real.buildMatrix(matrix.rowNum, matrix.colNum) { row, col ->
this - matrix[row, col] this@minus - matrix[row, col]
} }
// TODO: does this operation make sense? Should it be 'this/matrix[row, col]'? // TODO: does this operation make sense? Should it be 'this/matrix[row, col]'?
@ -87,29 +92,29 @@ public operator fun Double.minus(matrix: RealMatrix): RealMatrix =
@UnstableKMathAPI @UnstableKMathAPI
public operator fun RealMatrix.times(other: RealMatrix): RealMatrix = public operator fun RealMatrix.times(other: RealMatrix): RealMatrix =
MatrixContext.real.produce(rowNum, colNum) { row, col -> this[row, col] * other[row, col] } LinearSpace.real.buildMatrix(rowNum, colNum) { row, col -> this@times[row, col] * other[row, col] }
public operator fun RealMatrix.plus(other: RealMatrix): RealMatrix = public operator fun RealMatrix.plus(other: RealMatrix): RealMatrix =
MatrixContext.real.add(this, other) LinearSpace.real.run { this@plus + other }
public operator fun RealMatrix.minus(other: RealMatrix): RealMatrix = public operator fun RealMatrix.minus(other: RealMatrix): RealMatrix =
MatrixContext.real.produce(rowNum, colNum) { row, col -> this[row, col] - other[row, col] } LinearSpace.real.buildMatrix(rowNum, colNum) { row, col -> this@minus[row, col] - other[row, col] }
/* /*
* Operations on columns * Operations on columns
*/ */
public inline fun RealMatrix.appendColumn(crossinline mapper: (Buffer<Double>) -> Double): RealMatrix = public inline fun RealMatrix.appendColumn(crossinline mapper: (Buffer<Double>) -> Double): RealMatrix =
MatrixContext.real.produce(rowNum, colNum + 1) { row, col -> LinearSpace.real.buildMatrix(rowNum, colNum + 1) { row, col ->
if (col < colNum) if (col < colNum)
this[row, col] get(row, col)
else else
mapper(rows[row]) mapper(rows[row])
} }
public fun RealMatrix.extractColumns(columnRange: IntRange): RealMatrix = public fun RealMatrix.extractColumns(columnRange: IntRange): RealMatrix =
MatrixContext.real.produce(rowNum, columnRange.count()) { row, col -> LinearSpace.real.buildMatrix(rowNum, columnRange.count()) { row, col ->
this[row, columnRange.first + col] this@extractColumns[row, columnRange.first + col]
} }
public fun RealMatrix.extractColumn(columnIndex: Int): RealMatrix = public fun RealMatrix.extractColumn(columnIndex: Int): RealMatrix =
@ -141,14 +146,14 @@ public fun RealMatrix.max(): Double? = elements().map { (_, value) -> value }.ma
public fun RealMatrix.average(): Double = elements().map { (_, value) -> value }.average() public fun RealMatrix.average(): Double = elements().map { (_, value) -> value }.average()
public inline fun RealMatrix.map(crossinline transform: (Double) -> Double): RealMatrix = public inline fun RealMatrix.map(crossinline transform: (Double) -> Double): RealMatrix =
MatrixContext.real.produce(rowNum, colNum) { i, j -> LinearSpace.real.buildMatrix(rowNum, colNum) { i, j ->
transform(get(i, j)) transform(get(i, j))
} }
/** /**
* Inverse a square real matrix using LUP decomposition * Inverse a square real matrix using LUP decomposition
*/ */
public fun RealMatrix.inverseWithLup(): RealMatrix = MatrixContext.real.inverseWithLup(this) public fun RealMatrix.inverseWithLup(): RealMatrix = LinearSpace.real.inverseWithLup(this)
//extended operations //extended operations

@ -5,6 +5,7 @@ import space.kscience.kmath.operations.Norm
import space.kscience.kmath.structures.Buffer import space.kscience.kmath.structures.Buffer
import space.kscience.kmath.structures.asBuffer import space.kscience.kmath.structures.asBuffer
import space.kscience.kmath.structures.asIterable import space.kscience.kmath.structures.asIterable
import space.kscience.kmath.structures.indices
import kotlin.math.pow import kotlin.math.pow
import kotlin.math.sqrt import kotlin.math.sqrt
@ -30,7 +31,7 @@ public inline fun RealVector.mapIndexed(transform: (index: Int, value: Double) -
Buffer.real(size) { transform(it, get(it)) } Buffer.real(size) { transform(it, get(it)) }
public operator fun RealVector.plus(other: RealVector): RealVector { public operator fun RealVector.plus(other: RealVector): RealVector {
require(size == other.size){"Vector size $size expected but ${other.size} found"} require(size == other.size) { "Vector size $size expected but ${other.size} found" }
return mapIndexed { index, value -> value + other[index] } return mapIndexed { index, value -> value + other[index] }
} }
@ -41,7 +42,7 @@ public operator fun Number.plus(vector: RealVector): RealVector = vector + this
public operator fun RealVector.unaryMinus(): Buffer<Double> = map { -it } public operator fun RealVector.unaryMinus(): Buffer<Double> = map { -it }
public operator fun RealVector.minus(other: RealVector): RealVector { public operator fun RealVector.minus(other: RealVector): RealVector {
require(size == other.size){"Vector size $size expected but ${other.size} found"} require(size == other.size) { "Vector size $size expected but ${other.size} found" }
return mapIndexed { index, value -> value - other[index] } return mapIndexed { index, value -> value - other[index] }
} }
@ -50,7 +51,7 @@ public operator fun RealVector.minus(number: Number): RealVector = map { it - nu
public operator fun Number.minus(vector: RealVector): RealVector = vector.map { toDouble() - it } public operator fun Number.minus(vector: RealVector): RealVector = vector.map { toDouble() - it }
public operator fun RealVector.times(other: RealVector): RealVector { public operator fun RealVector.times(other: RealVector): RealVector {
require(size == other.size){"Vector size $size expected but ${other.size} found"} require(size == other.size) { "Vector size $size expected but ${other.size} found" }
return mapIndexed { index, value -> value * other[index] } return mapIndexed { index, value -> value * other[index] }
} }
@ -59,7 +60,7 @@ public operator fun RealVector.times(number: Number): RealVector = map { it * nu
public operator fun Number.times(vector: RealVector): RealVector = vector * this public operator fun Number.times(vector: RealVector): RealVector = vector * this
public operator fun RealVector.div(other: RealVector): RealVector { public operator fun RealVector.div(other: RealVector): RealVector {
require(size == other.size){"Vector size $size expected but ${other.size} found"} require(size == other.size) { "Vector size $size expected but ${other.size} found" }
return mapIndexed { index, value -> value / other[index] } return mapIndexed { index, value -> value / other[index] }
} }
@ -88,3 +89,13 @@ public fun tan(vector: RealVector): RealVector = vector.map { kotlin.math.tan(it
public fun ln(vector: RealVector): RealVector = vector.map { kotlin.math.ln(it) } public fun ln(vector: RealVector): RealVector = vector.map { kotlin.math.ln(it) }
public fun log10(vector: RealVector): RealVector = vector.map { kotlin.math.log10(it) } public fun log10(vector: RealVector): RealVector = vector.map { kotlin.math.log10(it) }
// reductions methods
public fun RealVector.sum(): Double {
var res = 0.0
for (i in indices) {
res += get(i)
}
return res
}

@ -1,31 +1,12 @@
package space.kscience.kmath.real package space.kscience.kmath.real
import space.kscience.kmath.linear.BufferMatrix import space.kscience.kmath.linear.LinearSpace
import space.kscience.kmath.structures.Buffer import space.kscience.kmath.linear.Matrix
import space.kscience.kmath.structures.RealBuffer
/** /**
* Optimized dot product for real matrices * Optimized dot product for real matrices
*/ */
public infix fun BufferMatrix<Double>.dot(other: BufferMatrix<Double>): BufferMatrix<Double> { public infix fun Matrix<Double>.dot(other: Matrix<Double>): Matrix<Double> = LinearSpace.real.run{
require(colNum == other.rowNum) { "Matrix dot operation dimension mismatch: ($rowNum, $colNum) x (${other.rowNum}, ${other.colNum})" } this@dot dot other
val resultArray = DoubleArray(this.rowNum * other.colNum)
//convert to array to insure there is no memory indirection
fun Buffer<Double>.unsafeArray() = if (this is RealBuffer)
this.array
else
DoubleArray(size) { get(it) }
val a = this.buffer.unsafeArray()
val b = other.buffer.unsafeArray()
for (i in (0 until rowNum))
for (j in (0 until other.colNum))
for (k in (0 until colNum))
resultArray[i * other.colNum + j] += a[i * colNum + k] * b[k * other.colNum + j]
val buffer = RealBuffer(resultArray)
return BufferMatrix(rowNum, other.colNum, buffer)
} }

@ -1,7 +1,7 @@
package kaceince.kmath.real package kaceince.kmath.real
import space.kscience.kmath.linear.Matrix import space.kscience.kmath.linear.LinearSpace
import space.kscience.kmath.linear.build import space.kscience.kmath.linear.matrix
import space.kscience.kmath.misc.UnstableKMathAPI import space.kscience.kmath.misc.UnstableKMathAPI
import space.kscience.kmath.real.* import space.kscience.kmath.real.*
import space.kscience.kmath.structures.contentEquals import space.kscience.kmath.structures.contentEquals
@ -9,6 +9,7 @@ import kotlin.test.Test
import kotlin.test.assertEquals import kotlin.test.assertEquals
import kotlin.test.assertTrue import kotlin.test.assertTrue
@UnstableKMathAPI
internal class RealMatrixTest { internal class RealMatrixTest {
@Test @Test
fun testSum() { fun testSum() {
@ -30,11 +31,11 @@ internal class RealMatrixTest {
@Test @Test
fun testRepeatStackVertical() { fun testRepeatStackVertical() {
val matrix1 = Matrix.build(2, 3)( val matrix1 = realMatrix(2, 3)(
1.0, 0.0, 0.0, 1.0, 0.0, 0.0,
0.0, 1.0, 2.0 0.0, 1.0, 2.0
) )
val matrix2 = Matrix.build(6, 3)( val matrix2 = realMatrix(6, 3)(
1.0, 0.0, 0.0, 1.0, 0.0, 0.0,
0.0, 1.0, 2.0, 0.0, 1.0, 2.0,
1.0, 0.0, 0.0, 1.0, 0.0, 0.0,
@ -47,12 +48,12 @@ internal class RealMatrixTest {
@Test @Test
fun testMatrixAndDouble() { fun testMatrixAndDouble() {
val matrix1 = Matrix.build(2, 3)( val matrix1 = realMatrix(2, 3)(
1.0, 0.0, 3.0, 1.0, 0.0, 3.0,
4.0, 6.0, 2.0 4.0, 6.0, 2.0
) )
val matrix2 = (matrix1 * 2.5 + 1.0 - 2.0) / 2.0 val matrix2 = (matrix1 * 2.5 + 1.0 - 2.0) / 2.0
val expectedResult = Matrix.build(2, 3)( val expectedResult = LinearSpace.real.matrix(2, 3)(
0.75, -0.5, 3.25, 0.75, -0.5, 3.25,
4.5, 7.0, 2.0 4.5, 7.0, 2.0
) )
@ -61,13 +62,13 @@ internal class RealMatrixTest {
@Test @Test
fun testDoubleAndMatrix() { fun testDoubleAndMatrix() {
val matrix1 = Matrix.build(2, 3)( val matrix1 = realMatrix(2, 3)(
1.0, 0.0, 3.0, 1.0, 0.0, 3.0,
4.0, 6.0, 2.0 4.0, 6.0, 2.0
) )
val matrix2 = 20.0 - (10.0 + (5.0 * matrix1)) val matrix2 = 20.0 - (10.0 + (5.0 * matrix1))
//val matrix2 = 10.0 + (5.0 * matrix1) //val matrix2 = 10.0 + (5.0 * matrix1)
val expectedResult = Matrix.build(2, 3)( val expectedResult = realMatrix(2, 3)(
5.0, 10.0, -5.0, 5.0, 10.0, -5.0,
-10.0, -20.0, 0.0 -10.0, -20.0, 0.0
) )
@ -76,15 +77,15 @@ internal class RealMatrixTest {
@Test @Test
fun testSquareAndPower() { fun testSquareAndPower() {
val matrix1 = Matrix.build(2, 3)( val matrix1 = realMatrix(2, 3)(
-1.0, 0.0, 3.0, -1.0, 0.0, 3.0,
4.0, -6.0, -2.0 4.0, -6.0, -2.0
) )
val matrix2 = Matrix.build(2, 3)( val matrix2 = realMatrix(2, 3)(
1.0, 0.0, 9.0, 1.0, 0.0, 9.0,
16.0, 36.0, 4.0 16.0, 36.0, 4.0
) )
val matrix3 = Matrix.build(2, 3)( val matrix3 = realMatrix(2, 3)(
-1.0, 0.0, 27.0, -1.0, 0.0, 27.0,
64.0, -216.0, -8.0 64.0, -216.0, -8.0
) )
@ -95,16 +96,16 @@ internal class RealMatrixTest {
@OptIn(UnstableKMathAPI::class) @OptIn(UnstableKMathAPI::class)
@Test @Test
fun testTwoMatrixOperations() { fun testTwoMatrixOperations() {
val matrix1 = Matrix.build(2, 3)( val matrix1 = realMatrix(2, 3)(
-1.0, 0.0, 3.0, -1.0, 0.0, 3.0,
4.0, -6.0, 7.0 4.0, -6.0, 7.0
) )
val matrix2 = Matrix.build(2, 3)( val matrix2 = realMatrix(2, 3)(
1.0, 0.0, 3.0, 1.0, 0.0, 3.0,
4.0, 6.0, -2.0 4.0, 6.0, -2.0
) )
val result = matrix1 * matrix2 + matrix1 - matrix2 val result = matrix1 * matrix2 + matrix1 - matrix2
val expectedResult = Matrix.build(2, 3)( val expectedResult = realMatrix(2, 3)(
-3.0, 0.0, 9.0, -3.0, 0.0, 9.0,
16.0, -48.0, -5.0 16.0, -48.0, -5.0
) )
@ -113,16 +114,16 @@ internal class RealMatrixTest {
@Test @Test
fun testColumnOperations() { fun testColumnOperations() {
val matrix1 = Matrix.build(2, 4)( val matrix1 = realMatrix(2, 4)(
-1.0, 0.0, 3.0, 15.0, -1.0, 0.0, 3.0, 15.0,
4.0, -6.0, 7.0, -11.0 4.0, -6.0, 7.0, -11.0
) )
val matrix2 = Matrix.build(2, 5)( val matrix2 = realMatrix(2, 5)(
-1.0, 0.0, 3.0, 15.0, -1.0, -1.0, 0.0, 3.0, 15.0, -1.0,
4.0, -6.0, 7.0, -11.0, 4.0 4.0, -6.0, 7.0, -11.0, 4.0
) )
val col1 = Matrix.build(2, 1)(0.0, -6.0) val col1 = realMatrix(2, 1)(0.0, -6.0)
val cols1to2 = Matrix.build(2, 2)( val cols1to2 = realMatrix(2, 2)(
0.0, 3.0, 0.0, 3.0,
-6.0, 7.0 -6.0, 7.0
) )
@ -147,7 +148,7 @@ internal class RealMatrixTest {
@Test @Test
fun testAllElementOperations() { fun testAllElementOperations() {
val matrix1 = Matrix.build(2, 4)( val matrix1 = LinearSpace.real.matrix(2, 4)(
-1.0, 0.0, 3.0, 15.0, -1.0, 0.0, 3.0, 15.0,
4.0, -6.0, 7.0, -11.0 4.0, -6.0, 7.0, -11.0
) )

@ -1,10 +1,8 @@
package kaceince.kmath.real package kaceince.kmath.real
import space.kscience.kmath.linear.MatrixContext import space.kscience.kmath.linear.LinearSpace
import space.kscience.kmath.linear.asMatrix import space.kscience.kmath.linear.asMatrix
import space.kscience.kmath.linear.real
import space.kscience.kmath.linear.transpose import space.kscience.kmath.linear.transpose
import space.kscience.kmath.operations.invoke
import space.kscience.kmath.real.plus import space.kscience.kmath.real.plus
import space.kscience.kmath.structures.Buffer import space.kscience.kmath.structures.Buffer
import kotlin.test.Test import kotlin.test.Test
@ -32,7 +30,7 @@ internal class RealVectorTest {
val vector2 = Buffer.real(5) { 5 - it.toDouble() } val vector2 = Buffer.real(5) { 5 - it.toDouble() }
val matrix1 = vector1.asMatrix() val matrix1 = vector1.asMatrix()
val matrix2 = vector2.asMatrix().transpose() val matrix2 = vector2.asMatrix().transpose()
val product = MatrixContext.real { matrix1 dot matrix2 } val product = LinearSpace.real.run { matrix1 dot matrix2 }
assertEquals(5.0, product[1, 0]) assertEquals(5.0, product[1, 0])
assertEquals(6.0, product[2, 2]) assertEquals(6.0, product[2, 2])
} }

@ -1,7 +1,8 @@
package space.kscience.kmath.functions package space.kscience.kmath.functions
import space.kscience.kmath.operations.Group
import space.kscience.kmath.operations.Ring import space.kscience.kmath.operations.Ring
import space.kscience.kmath.operations.Space import space.kscience.kmath.operations.ScaleOperations
import space.kscience.kmath.operations.invoke import space.kscience.kmath.operations.invoke
import kotlin.contracts.InvocationKind import kotlin.contracts.InvocationKind
import kotlin.contracts.contract import kotlin.contracts.contract
@ -41,9 +42,15 @@ public fun <T : Any, C : Ring<T>> Polynomial<T>.asFunction(ring: C): (T) -> T =
/** /**
* An algebra for polynomials * An algebra for polynomials
*/ */
public class PolynomialSpace<T : Any, C : Ring<T>>(private val ring: C) : Space<Polynomial<T>> { public class PolynomialSpace<T : Any, C>(
private val ring: C,
) : Group<Polynomial<T>>, ScaleOperations<Polynomial<T>> where C : Ring<T>, C : ScaleOperations<T> {
public override val zero: Polynomial<T> = Polynomial(emptyList()) public override val zero: Polynomial<T> = Polynomial(emptyList())
override fun Polynomial<T>.unaryMinus(): Polynomial<T> = with(ring) {
Polynomial(coefficients.map { -it })
}
public override fun add(a: Polynomial<T>, b: Polynomial<T>): Polynomial<T> { public override fun add(a: Polynomial<T>, b: Polynomial<T>): Polynomial<T> {
val dim = max(a.coefficients.size, b.coefficients.size) val dim = max(a.coefficients.size, b.coefficients.size)
@ -54,13 +61,13 @@ public class PolynomialSpace<T : Any, C : Ring<T>>(private val ring: C) : Space<
} }
} }
public override fun multiply(a: Polynomial<T>, k: Number): Polynomial<T> = public override fun scale(a: Polynomial<T>, value: Double): Polynomial<T> =
ring { Polynomial(List(a.coefficients.size) { index -> a.coefficients[index] * k }) } ring { Polynomial(List(a.coefficients.size) { index -> a.coefficients[index] * value }) }
public operator fun Polynomial<T>.invoke(arg: T): T = value(ring, arg) public operator fun Polynomial<T>.invoke(arg: T): T = value(ring, arg)
} }
public inline fun <T : Any, C : Ring<T>, R> C.polynomial(block: PolynomialSpace<T, C>.() -> R): R { public inline fun <T : Any, C, R> C.polynomial(block: PolynomialSpace<T, C>.() -> R): R where C : Ring<T>, C : ScaleOperations<T> {
contract { callsInPlace(block, InvocationKind.EXACTLY_ONCE) } contract { callsInPlace(block, InvocationKind.EXACTLY_ONCE) }
return PolynomialSpace(this).block() return PolynomialSpace(this).block()
} }

@ -13,7 +13,7 @@ import space.kscience.kmath.structures.MutableBufferFactory
*/ */
public class SplineInterpolator<T : Comparable<T>>( public class SplineInterpolator<T : Comparable<T>>(
public override val algebra: Field<T>, public override val algebra: Field<T>,
public val bufferFactory: MutableBufferFactory<T> public val bufferFactory: MutableBufferFactory<T>,
) : PolynomialInterpolator<T> { ) : PolynomialInterpolator<T> {
//TODO possibly optimize zeroed buffers //TODO possibly optimize zeroed buffers

Some files were not shown because too many files have changed in this diff Show More