refactoring #303
@ -53,14 +53,15 @@ public object DoubleLinearOpsTensorAlgebra :
|
|||||||
val lTensor = luTensor.zeroesLike()
|
val lTensor = luTensor.zeroesLike()
|
||||||
val uTensor = luTensor.zeroesLike()
|
val uTensor = luTensor.zeroesLike()
|
||||||
|
|
||||||
for ((pairLU, lu) in lTensor.matrixSequence().zip(uTensor.matrixSequence())
|
lTensor.matrixSequence()
|
||||||
.zip(luTensor.tensor.matrixSequence())) {
|
.zip(uTensor.matrixSequence())
|
||||||
|
.zip(luTensor.tensor.matrixSequence())
|
||||||
|
.forEach { (pairLU, lu) ->
|
||||||
val (l, u) = pairLU
|
val (l, u) = pairLU
|
||||||
luPivotHelper(l.as2D(), u.as2D(), lu.as2D(), n)
|
luPivotHelper(l.as2D(), u.as2D(), lu.as2D(), n)
|
||||||
}
|
}
|
||||||
|
|
||||||
return Triple(pTensor, lTensor, uTensor)
|
return Triple(pTensor, lTensor, uTensor)
|
||||||
|
|
||||||
}
|
}
|
||||||
|
|
||||||
public fun TensorStructure<Double>.cholesky(epsilon: Double): DoubleTensor {
|
public fun TensorStructure<Double>.cholesky(epsilon: Double): DoubleTensor {
|
||||||
@ -82,11 +83,13 @@ public object DoubleLinearOpsTensorAlgebra :
|
|||||||
checkSquareMatrix(shape)
|
checkSquareMatrix(shape)
|
||||||
val qTensor = zeroesLike()
|
val qTensor = zeroesLike()
|
||||||
val rTensor = zeroesLike()
|
val rTensor = zeroesLike()
|
||||||
val seq = tensor.matrixSequence().zip((qTensor.matrixSequence().zip(rTensor.matrixSequence())))
|
tensor.matrixSequence()
|
||||||
for ((matrix, qr) in seq) {
|
.zip((qTensor.matrixSequence()
|
||||||
|
.zip(rTensor.matrixSequence()))).forEach { (matrix, qr) ->
|
||||||
val (q, r) = qr
|
val (q, r) = qr
|
||||||
qrHelper(matrix.asTensor(), q.asTensor(), r.as2D())
|
qrHelper(matrix.asTensor(), q.asTensor(), r.as2D())
|
||||||
}
|
}
|
||||||
|
|
||||||
return qTensor to rTensor
|
return qTensor to rTensor
|
||||||
}
|
}
|
||||||
|
|
||||||
@ -97,20 +100,24 @@ public object DoubleLinearOpsTensorAlgebra :
|
|||||||
val size = tensor.linearStructure.dim
|
val size = tensor.linearStructure.dim
|
||||||
val commonShape = tensor.shape.sliceArray(0 until size - 2)
|
val commonShape = tensor.shape.sliceArray(0 until size - 2)
|
||||||
val (n, m) = tensor.shape.sliceArray(size - 2 until size)
|
val (n, m) = tensor.shape.sliceArray(size - 2 until size)
|
||||||
val resU = zeros(commonShape + intArrayOf(min(n, m), n))
|
val uTensor = zeros(commonShape + intArrayOf(min(n, m), n))
|
||||||
val resS = zeros(commonShape + intArrayOf(min(n, m)))
|
val sTensor = zeros(commonShape + intArrayOf(min(n, m)))
|
||||||
val resV = zeros(commonShape + intArrayOf(min(n, m), m))
|
val vTensor = zeros(commonShape + intArrayOf(min(n, m), m))
|
||||||
|
|
||||||
for ((matrix, USV) in tensor.matrixSequence()
|
tensor.matrixSequence()
|
||||||
.zip(resU.matrixSequence().zip(resS.vectorSequence().zip(resV.matrixSequence())))) {
|
.zip(uTensor.matrixSequence()
|
||||||
|
.zip(sTensor.vectorSequence()
|
||||||
|
.zip(vTensor.matrixSequence()))).forEach { (matrix, USV) ->
|
||||||
val matrixSize = matrix.shape.reduce { acc, i -> acc * i }
|
val matrixSize = matrix.shape.reduce { acc, i -> acc * i }
|
||||||
val curMatrix = DoubleTensor(
|
val curMatrix = DoubleTensor(
|
||||||
matrix.shape,
|
matrix.shape,
|
||||||
matrix.mutableBuffer.array().slice(matrix.bufferStart until matrix.bufferStart + matrixSize).toDoubleArray()
|
matrix.mutableBuffer.array().slice(matrix.bufferStart until matrix.bufferStart + matrixSize)
|
||||||
|
.toDoubleArray()
|
||||||
)
|
)
|
||||||
svdHelper(curMatrix, USV, m, n, epsilon)
|
svdHelper(curMatrix, USV, m, n, epsilon)
|
||||||
}
|
}
|
||||||
return Triple(resU.transpose(), resS, resV.transpose())
|
|
||||||
|
return Triple(uTensor.transpose(), sTensor, vTensor.transpose())
|
||||||
}
|
}
|
||||||
|
|
||||||
override fun TensorStructure<Double>.symEig(): Pair<DoubleTensor, DoubleTensor> =
|
override fun TensorStructure<Double>.symEig(): Pair<DoubleTensor, DoubleTensor> =
|
||||||
@ -127,7 +134,7 @@ public object DoubleLinearOpsTensorAlgebra :
|
|||||||
cleanSymHelper(matrix.as2D(), n)
|
cleanSymHelper(matrix.as2D(), n)
|
||||||
|
|
||||||
val eig = (utv dot s.view(shp)).view(s.shape)
|
val eig = (utv dot s.view(shp)).view(s.shape)
|
||||||
return Pair(eig, v)
|
return eig to v
|
||||||
}
|
}
|
||||||
|
|
||||||
public fun TensorStructure<Double>.detLU(epsilon: Double = 1e-9): DoubleTensor {
|
public fun TensorStructure<Double>.detLU(epsilon: Double = 1e-9): DoubleTensor {
|
||||||
|
@ -2,7 +2,7 @@ package space.kscience.kmath.tensors.core
|
|||||||
|
|
||||||
import kotlin.math.max
|
import kotlin.math.max
|
||||||
|
|
||||||
internal inline fun multiIndexBroadCasting(tensor: DoubleTensor, resTensor: DoubleTensor, linearSize: Int) {
|
internal fun multiIndexBroadCasting(tensor: DoubleTensor, resTensor: DoubleTensor, linearSize: Int) {
|
||||||
for (linearIndex in 0 until linearSize) {
|
for (linearIndex in 0 until linearSize) {
|
||||||
val totalMultiIndex = resTensor.linearStructure.index(linearIndex)
|
val totalMultiIndex = resTensor.linearStructure.index(linearIndex)
|
||||||
val curMultiIndex = tensor.shape.copyOf()
|
val curMultiIndex = tensor.shape.copyOf()
|
||||||
@ -23,7 +23,7 @@ internal inline fun multiIndexBroadCasting(tensor: DoubleTensor, resTensor: Doub
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun broadcastShapes(vararg shapes: IntArray): IntArray {
|
internal fun broadcastShapes(vararg shapes: IntArray): IntArray {
|
||||||
var totalDim = 0
|
var totalDim = 0
|
||||||
for (shape in shapes) {
|
for (shape in shapes) {
|
||||||
totalDim = max(totalDim, shape.size)
|
totalDim = max(totalDim, shape.size)
|
||||||
@ -51,7 +51,7 @@ internal inline fun broadcastShapes(vararg shapes: IntArray): IntArray {
|
|||||||
return totalShape
|
return totalShape
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun broadcastTo(tensor: DoubleTensor, newShape: IntArray): DoubleTensor {
|
internal fun broadcastTo(tensor: DoubleTensor, newShape: IntArray): DoubleTensor {
|
||||||
if (tensor.shape.size > newShape.size) {
|
if (tensor.shape.size > newShape.size) {
|
||||||
throw RuntimeException("Tensor is not compatible with the new shape")
|
throw RuntimeException("Tensor is not compatible with the new shape")
|
||||||
}
|
}
|
||||||
@ -71,7 +71,7 @@ internal inline fun broadcastTo(tensor: DoubleTensor, newShape: IntArray): Doubl
|
|||||||
return resTensor
|
return resTensor
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun broadcastTensors(vararg tensors: DoubleTensor): List<DoubleTensor> {
|
internal fun broadcastTensors(vararg tensors: DoubleTensor): List<DoubleTensor> {
|
||||||
val totalShape = broadcastShapes(*(tensors.map { it.shape }).toTypedArray())
|
val totalShape = broadcastShapes(*(tensors.map { it.shape }).toTypedArray())
|
||||||
val n = totalShape.reduce { acc, i -> acc * i }
|
val n = totalShape.reduce { acc, i -> acc * i }
|
||||||
|
|
||||||
@ -85,7 +85,7 @@ internal inline fun broadcastTensors(vararg tensors: DoubleTensor): List<DoubleT
|
|||||||
return res
|
return res
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun broadcastOuterTensors(vararg tensors: DoubleTensor): List<DoubleTensor> {
|
internal fun broadcastOuterTensors(vararg tensors: DoubleTensor): List<DoubleTensor> {
|
||||||
val onlyTwoDims = tensors.asSequence().onEach {
|
val onlyTwoDims = tensors.asSequence().onEach {
|
||||||
require(it.shape.size >= 2) {
|
require(it.shape.size >= 2) {
|
||||||
throw RuntimeException("Tensors must have at least 2 dimensions")
|
throw RuntimeException("Tensors must have at least 2 dimensions")
|
||||||
@ -99,7 +99,7 @@ internal inline fun broadcastOuterTensors(vararg tensors: DoubleTensor): List<Do
|
|||||||
val totalShape = broadcastShapes(*(tensors.map { it.shape.sliceArray(0..it.shape.size - 3) }).toTypedArray())
|
val totalShape = broadcastShapes(*(tensors.map { it.shape.sliceArray(0..it.shape.size - 3) }).toTypedArray())
|
||||||
val n = totalShape.reduce { acc, i -> acc * i }
|
val n = totalShape.reduce { acc, i -> acc * i }
|
||||||
|
|
||||||
val res = ArrayList<DoubleTensor>(0)
|
return buildList {
|
||||||
for (tensor in tensors) {
|
for (tensor in tensors) {
|
||||||
val matrixShape = tensor.shape.sliceArray(tensor.shape.size - 2 until tensor.shape.size).copyOf()
|
val matrixShape = tensor.shape.sliceArray(tensor.shape.size - 2 until tensor.shape.size).copyOf()
|
||||||
val matrixSize = matrixShape[0] * matrixShape[1]
|
val matrixSize = matrixShape[0] * matrixShape[1]
|
||||||
@ -137,8 +137,7 @@ internal inline fun broadcastOuterTensors(vararg tensors: DoubleTensor): List<Do
|
|||||||
newTensor.mutableBuffer.array()[newTensor.bufferStart + curLinearIndex]
|
newTensor.mutableBuffer.array()[newTensor.bufferStart + curLinearIndex]
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
res += resTensor
|
add(resTensor)
|
||||||
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
return res
|
|
||||||
}
|
}
|
@ -5,38 +5,38 @@ import space.kscience.kmath.tensors.core.algebras.DoubleLinearOpsTensorAlgebra
|
|||||||
import space.kscience.kmath.tensors.core.algebras.DoubleTensorAlgebra
|
import space.kscience.kmath.tensors.core.algebras.DoubleTensorAlgebra
|
||||||
|
|
||||||
|
|
||||||
internal inline fun checkEmptyShape(shape: IntArray): Unit =
|
internal fun checkEmptyShape(shape: IntArray): Unit =
|
||||||
check(shape.isNotEmpty()) {
|
check(shape.isNotEmpty()) {
|
||||||
"Illegal empty shape provided"
|
"Illegal empty shape provided"
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun checkEmptyDoubleBuffer(buffer: DoubleArray): Unit =
|
internal fun checkEmptyDoubleBuffer(buffer: DoubleArray): Unit =
|
||||||
check(buffer.isNotEmpty()) {
|
check(buffer.isNotEmpty()) {
|
||||||
"Illegal empty buffer provided"
|
"Illegal empty buffer provided"
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun checkBufferShapeConsistency(shape: IntArray, buffer: DoubleArray): Unit =
|
internal fun checkBufferShapeConsistency(shape: IntArray, buffer: DoubleArray): Unit =
|
||||||
check(buffer.size == shape.reduce(Int::times)) {
|
check(buffer.size == shape.reduce(Int::times)) {
|
||||||
"Inconsistent shape ${shape.toList()} for buffer of size ${buffer.size} provided"
|
"Inconsistent shape ${shape.toList()} for buffer of size ${buffer.size} provided"
|
||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
internal inline fun <T> checkShapesCompatible(a: TensorStructure<T>, b: TensorStructure<T>): Unit =
|
internal fun <T> checkShapesCompatible(a: TensorStructure<T>, b: TensorStructure<T>): Unit =
|
||||||
check(a.shape contentEquals b.shape) {
|
check(a.shape contentEquals b.shape) {
|
||||||
"Incompatible shapes ${a.shape.toList()} and ${b.shape.toList()} "
|
"Incompatible shapes ${a.shape.toList()} and ${b.shape.toList()} "
|
||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
internal inline fun checkTranspose(dim: Int, i: Int, j: Int): Unit =
|
internal fun checkTranspose(dim: Int, i: Int, j: Int): Unit =
|
||||||
check((i < dim) and (j < dim)) {
|
check((i < dim) and (j < dim)) {
|
||||||
"Cannot transpose $i to $j for a tensor of dim $dim"
|
"Cannot transpose $i to $j for a tensor of dim $dim"
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun <T> checkView(a: TensorStructure<T>, shape: IntArray): Unit =
|
internal fun <T> checkView(a: TensorStructure<T>, shape: IntArray): Unit =
|
||||||
check(a.shape.reduce(Int::times) == shape.reduce(Int::times))
|
check(a.shape.reduce(Int::times) == shape.reduce(Int::times))
|
||||||
|
|
||||||
|
|
||||||
internal inline fun checkSquareMatrix(shape: IntArray): Unit {
|
internal fun checkSquareMatrix(shape: IntArray): Unit {
|
||||||
val n = shape.size
|
val n = shape.size
|
||||||
check(n >= 2) {
|
check(n >= 2) {
|
||||||
"Expected tensor with 2 or more dimensions, got size $n instead"
|
"Expected tensor with 2 or more dimensions, got size $n instead"
|
||||||
@ -46,14 +46,14 @@ internal inline fun checkSquareMatrix(shape: IntArray): Unit {
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun DoubleTensorAlgebra.checkSymmetric(
|
internal fun DoubleTensorAlgebra.checkSymmetric(
|
||||||
tensor: TensorStructure<Double>, epsilon: Double = 1e-6
|
tensor: TensorStructure<Double>, epsilon: Double = 1e-6
|
||||||
): Unit =
|
): Unit =
|
||||||
check(tensor.eq(tensor.transpose(), epsilon)) {
|
check(tensor.eq(tensor.transpose(), epsilon)) {
|
||||||
"Tensor is not symmetric about the last 2 dimensions at precision $epsilon"
|
"Tensor is not symmetric about the last 2 dimensions at precision $epsilon"
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun DoubleLinearOpsTensorAlgebra.checkPositiveDefinite(
|
internal fun DoubleLinearOpsTensorAlgebra.checkPositiveDefinite(
|
||||||
tensor: DoubleTensor, epsilon: Double = 1e-6
|
tensor: DoubleTensor, epsilon: Double = 1e-6
|
||||||
): Unit {
|
): Unit {
|
||||||
checkSymmetric(tensor, epsilon)
|
checkSymmetric(tensor, epsilon)
|
||||||
|
@ -13,7 +13,7 @@ import kotlin.math.sign
|
|||||||
import kotlin.math.sqrt
|
import kotlin.math.sqrt
|
||||||
|
|
||||||
|
|
||||||
internal inline fun <T> BufferedTensor<T>.vectorSequence(): Sequence<BufferedTensor<T>> = sequence {
|
internal fun <T> BufferedTensor<T>.vectorSequence(): Sequence<BufferedTensor<T>> = sequence {
|
||||||
val n = shape.size
|
val n = shape.size
|
||||||
val vectorOffset = shape[n - 1]
|
val vectorOffset = shape[n - 1]
|
||||||
val vectorShape = intArrayOf(shape.last())
|
val vectorShape = intArrayOf(shape.last())
|
||||||
@ -23,9 +23,9 @@ internal inline fun <T> BufferedTensor<T>.vectorSequence(): Sequence<BufferedTen
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun <T> BufferedTensor<T>.matrixSequence(): Sequence<BufferedTensor<T>> = sequence {
|
internal fun <T> BufferedTensor<T>.matrixSequence(): Sequence<BufferedTensor<T>> = sequence {
|
||||||
check(shape.size >= 2) { "todo" }
|
|
||||||
val n = shape.size
|
val n = shape.size
|
||||||
|
check(n >= 2) { "Expected tensor with 2 or more dimensions, got size $n" }
|
||||||
val matrixOffset = shape[n - 1] * shape[n - 2]
|
val matrixOffset = shape[n - 1] * shape[n - 2]
|
||||||
val matrixShape = intArrayOf(shape[n - 2], shape[n - 1])
|
val matrixShape = intArrayOf(shape[n - 2], shape[n - 1])
|
||||||
for (offset in 0 until numElements step matrixOffset) {
|
for (offset in 0 until numElements step matrixOffset) {
|
||||||
@ -46,8 +46,7 @@ internal inline fun <T> BufferedTensor<T>.forEachMatrix(matrixAction: (BufferedT
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
|
internal fun dotHelper(
|
||||||
internal inline fun dotHelper(
|
|
||||||
a: MutableStructure2D<Double>,
|
a: MutableStructure2D<Double>,
|
||||||
b: MutableStructure2D<Double>,
|
b: MutableStructure2D<Double>,
|
||||||
res: MutableStructure2D<Double>,
|
res: MutableStructure2D<Double>,
|
||||||
@ -64,10 +63,11 @@ internal inline fun dotHelper(
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun luHelper(
|
internal fun luHelper(
|
||||||
lu: MutableStructure2D<Double>,
|
lu: MutableStructure2D<Double>,
|
||||||
pivots: MutableStructure1D<Int>,
|
pivots: MutableStructure1D<Int>,
|
||||||
epsilon: Double): Boolean {
|
epsilon: Double
|
||||||
|
): Boolean {
|
||||||
|
|
||||||
val m = lu.rowNum
|
val m = lu.rowNum
|
||||||
|
|
||||||
@ -114,7 +114,7 @@ internal inline fun luHelper(
|
|||||||
return false
|
return false
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun <T> BufferedTensor<T>.setUpPivots(): IntTensor {
|
internal fun <T> BufferedTensor<T>.setUpPivots(): IntTensor {
|
||||||
val n = this.shape.size
|
val n = this.shape.size
|
||||||
val m = this.shape.last()
|
val m = this.shape.last()
|
||||||
val pivotsShape = IntArray(n - 1) { i -> this.shape[i] }
|
val pivotsShape = IntArray(n - 1) { i -> this.shape[i] }
|
||||||
@ -126,9 +126,10 @@ internal inline fun <T> BufferedTensor<T>.setUpPivots(): IntTensor {
|
|||||||
)
|
)
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun DoubleLinearOpsTensorAlgebra.computeLU(
|
internal fun DoubleLinearOpsTensorAlgebra.computeLU(
|
||||||
tensor: DoubleTensor,
|
tensor: DoubleTensor,
|
||||||
epsilon: Double): Pair<DoubleTensor, IntTensor>? {
|
epsilon: Double
|
||||||
|
): Pair<DoubleTensor, IntTensor>? {
|
||||||
|
|
||||||
checkSquareMatrix(tensor.shape)
|
checkSquareMatrix(tensor.shape)
|
||||||
val luTensor = tensor.copy()
|
val luTensor = tensor.copy()
|
||||||
@ -141,7 +142,7 @@ internal inline fun DoubleLinearOpsTensorAlgebra.computeLU(
|
|||||||
return Pair(luTensor, pivotsTensor)
|
return Pair(luTensor, pivotsTensor)
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun pivInit(
|
internal fun pivInit(
|
||||||
p: MutableStructure2D<Double>,
|
p: MutableStructure2D<Double>,
|
||||||
pivot: MutableStructure1D<Int>,
|
pivot: MutableStructure1D<Int>,
|
||||||
n: Int
|
n: Int
|
||||||
@ -151,7 +152,7 @@ internal inline fun pivInit(
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun luPivotHelper(
|
internal fun luPivotHelper(
|
||||||
l: MutableStructure2D<Double>,
|
l: MutableStructure2D<Double>,
|
||||||
u: MutableStructure2D<Double>,
|
u: MutableStructure2D<Double>,
|
||||||
lu: MutableStructure2D<Double>,
|
lu: MutableStructure2D<Double>,
|
||||||
@ -172,7 +173,7 @@ internal inline fun luPivotHelper(
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun choleskyHelper(
|
internal fun choleskyHelper(
|
||||||
a: MutableStructure2D<Double>,
|
a: MutableStructure2D<Double>,
|
||||||
l: MutableStructure2D<Double>,
|
l: MutableStructure2D<Double>,
|
||||||
n: Int
|
n: Int
|
||||||
@ -193,7 +194,7 @@ internal inline fun choleskyHelper(
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun luMatrixDet(lu: MutableStructure2D<Double>, pivots: MutableStructure1D<Int>): Double {
|
internal fun luMatrixDet(lu: MutableStructure2D<Double>, pivots: MutableStructure1D<Int>): Double {
|
||||||
if (lu[0, 0] == 0.0) {
|
if (lu[0, 0] == 0.0) {
|
||||||
return 0.0
|
return 0.0
|
||||||
}
|
}
|
||||||
@ -202,7 +203,7 @@ internal inline fun luMatrixDet(lu: MutableStructure2D<Double>, pivots: MutableS
|
|||||||
return (0 until m).asSequence().map { lu[it, it] }.fold(sign) { left, right -> left * right }
|
return (0 until m).asSequence().map { lu[it, it] }.fold(sign) { left, right -> left * right }
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun luMatrixInv(
|
internal fun luMatrixInv(
|
||||||
lu: MutableStructure2D<Double>,
|
lu: MutableStructure2D<Double>,
|
||||||
pivots: MutableStructure1D<Int>,
|
pivots: MutableStructure1D<Int>,
|
||||||
invMatrix: MutableStructure2D<Double>
|
invMatrix: MutableStructure2D<Double>
|
||||||
@ -229,7 +230,7 @@ internal inline fun luMatrixInv(
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun DoubleLinearOpsTensorAlgebra.qrHelper(
|
internal fun DoubleLinearOpsTensorAlgebra.qrHelper(
|
||||||
matrix: DoubleTensor,
|
matrix: DoubleTensor,
|
||||||
q: DoubleTensor,
|
q: DoubleTensor,
|
||||||
r: MutableStructure2D<Double>
|
r: MutableStructure2D<Double>
|
||||||
@ -259,7 +260,7 @@ internal inline fun DoubleLinearOpsTensorAlgebra.qrHelper(
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun DoubleLinearOpsTensorAlgebra.svd1d(a: DoubleTensor, epsilon: Double = 1e-10): DoubleTensor {
|
internal fun DoubleLinearOpsTensorAlgebra.svd1d(a: DoubleTensor, epsilon: Double = 1e-10): DoubleTensor {
|
||||||
val (n, m) = a.shape
|
val (n, m) = a.shape
|
||||||
var v: DoubleTensor
|
var v: DoubleTensor
|
||||||
val b: DoubleTensor
|
val b: DoubleTensor
|
||||||
@ -283,7 +284,7 @@ internal inline fun DoubleLinearOpsTensorAlgebra.svd1d(a: DoubleTensor, epsilon:
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun DoubleLinearOpsTensorAlgebra.svdHelper(
|
internal fun DoubleLinearOpsTensorAlgebra.svdHelper(
|
||||||
matrix: DoubleTensor,
|
matrix: DoubleTensor,
|
||||||
USV: Pair<BufferedTensor<Double>, Pair<BufferedTensor<Double>, BufferedTensor<Double>>>,
|
USV: Pair<BufferedTensor<Double>, Pair<BufferedTensor<Double>, BufferedTensor<Double>>>,
|
||||||
m: Int, n: Int, epsilon: Double
|
m: Int, n: Int, epsilon: Double
|
||||||
@ -335,7 +336,7 @@ internal inline fun DoubleLinearOpsTensorAlgebra.svdHelper(
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun cleanSymHelper(matrix: MutableStructure2D<Double>, n: Int) {
|
internal fun cleanSymHelper(matrix: MutableStructure2D<Double>, n: Int) {
|
||||||
for (i in 0 until n)
|
for (i in 0 until n)
|
||||||
for (j in 0 until n) {
|
for (j in 0 until n) {
|
||||||
if (i == j) {
|
if (i == j) {
|
||||||
|
@ -23,19 +23,19 @@ internal fun Buffer<Double>.array(): DoubleArray = when (this) {
|
|||||||
else -> this.toDoubleArray()
|
else -> this.toDoubleArray()
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun getRandomNormals(n: Int, seed: Long): DoubleArray {
|
internal fun getRandomNormals(n: Int, seed: Long): DoubleArray {
|
||||||
val distribution = GaussianSampler(0.0, 1.0)
|
val distribution = GaussianSampler(0.0, 1.0)
|
||||||
val generator = RandomGenerator.default(seed)
|
val generator = RandomGenerator.default(seed)
|
||||||
return distribution.sample(generator).nextBufferBlocking(n).toDoubleArray()
|
return distribution.sample(generator).nextBufferBlocking(n).toDoubleArray()
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun getRandomUnitVector(n: Int, seed: Long): DoubleArray {
|
internal fun getRandomUnitVector(n: Int, seed: Long): DoubleArray {
|
||||||
val unnorm = getRandomNormals(n, seed)
|
val unnorm = getRandomNormals(n, seed)
|
||||||
val norm = sqrt(unnorm.map { it * it }.sum())
|
val norm = sqrt(unnorm.map { it * it }.sum())
|
||||||
return unnorm.map { it / norm }.toDoubleArray()
|
return unnorm.map { it / norm }.toDoubleArray()
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun minusIndexFrom(n: Int, i: Int) : Int = if (i >= 0) i else {
|
internal fun minusIndexFrom(n: Int, i: Int): Int = if (i >= 0) i else {
|
||||||
val ii = n + i
|
val ii = n + i
|
||||||
check(ii >= 0) {
|
check(ii >= 0) {
|
||||||
"Out of bound index $i for tensor of dim $n"
|
"Out of bound index $i for tensor of dim $n"
|
||||||
@ -43,15 +43,16 @@ internal inline fun minusIndexFrom(n: Int, i: Int) : Int = if (i >= 0) i else {
|
|||||||
ii
|
ii
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun <T> BufferedTensor<T>.minusIndex(i: Int): Int = minusIndexFrom(this.dimension, i)
|
internal fun <T> BufferedTensor<T>.minusIndex(i: Int): Int = minusIndexFrom(this.dimension, i)
|
||||||
|
|
||||||
internal inline fun format(value: Double, digits: Int = 4): String {
|
internal fun format(value: Double, digits: Int = 4): String {
|
||||||
val ten = 10.0
|
val ten = 10.0
|
||||||
val approxOrder = if (value == 0.0) 0 else ceil(log10(abs(value))).toInt()
|
val approxOrder = if (value == 0.0) 0 else ceil(log10(abs(value))).toInt()
|
||||||
val order = if (
|
val order = if (
|
||||||
((value % ten) == 0.0) or
|
((value % ten) == 0.0) or
|
||||||
(value == 1.0) or
|
(value == 1.0) or
|
||||||
((1/value) % ten == 0.0)) approxOrder else approxOrder - 1
|
((1 / value) % ten == 0.0)
|
||||||
|
) approxOrder else approxOrder - 1
|
||||||
val lead = value / ten.pow(order)
|
val lead = value / ten.pow(order)
|
||||||
val leadDisplay = round(lead * ten.pow(digits)) / ten.pow(digits)
|
val leadDisplay = round(lead * ten.pow(digits)) / ten.pow(digits)
|
||||||
val orderDisplay = if (order == 0) "" else if (order > 0) "E+$order" else "E$order"
|
val orderDisplay = if (order == 0) "" else if (order > 0) "E+$order" else "E$order"
|
||||||
@ -63,7 +64,7 @@ internal inline fun format(value: Double, digits: Int = 4): String {
|
|||||||
return "$res$endSpace"
|
return "$res$endSpace"
|
||||||
}
|
}
|
||||||
|
|
||||||
internal inline fun DoubleTensor.toPrettyString(): String = buildString {
|
internal fun DoubleTensor.toPrettyString(): String = buildString {
|
||||||
var offset = 0
|
var offset = 0
|
||||||
val shape = this@toPrettyString.shape
|
val shape = this@toPrettyString.shape
|
||||||
val linearStructure = this@toPrettyString.linearStructure
|
val linearStructure = this@toPrettyString.linearStructure
|
||||||
@ -72,32 +73,36 @@ internal inline fun DoubleTensor.toPrettyString(): String = buildString {
|
|||||||
append(initString)
|
append(initString)
|
||||||
var charOffset = 3
|
var charOffset = 3
|
||||||
for (vector in vectorSequence()) {
|
for (vector in vectorSequence()) {
|
||||||
append(" ".repeat(charOffset))
|
repeat(charOffset) { append(' ') }
|
||||||
val index = linearStructure.index(offset)
|
val index = linearStructure.index(offset)
|
||||||
for (ind in index.reversed()) {
|
for (ind in index.reversed()) {
|
||||||
if (ind != 0) {
|
if (ind != 0) {
|
||||||
break
|
break
|
||||||
}
|
}
|
||||||
append("[")
|
append('[')
|
||||||
charOffset += 1
|
charOffset += 1
|
||||||
}
|
}
|
||||||
|
|
||||||
val values = vector.as1D().toMutableList().map(::format)
|
val values = vector.as1D().toMutableList().map(::format)
|
||||||
|
|
||||||
append(values.joinToString(", "))
|
values.joinTo(this, separator = ", ")
|
||||||
append("]")
|
|
||||||
|
append(']')
|
||||||
charOffset -= 1
|
charOffset -= 1
|
||||||
for ((ind, maxInd) in index.reversed().zip(shape.reversed()).drop(1)){
|
|
||||||
|
index.reversed().zip(shape.reversed()).drop(1).forEach { (ind, maxInd) ->
|
||||||
if (ind != maxInd - 1) {
|
if (ind != maxInd - 1) {
|
||||||
break
|
return@forEach
|
||||||
}
|
}
|
||||||
append("]")
|
append(']')
|
||||||
charOffset -= 1
|
charOffset -= 1
|
||||||
}
|
}
|
||||||
|
|
||||||
offset += vectorSize
|
offset += vectorSize
|
||||||
if (this@toPrettyString.numElements == offset) {
|
if (this@toPrettyString.numElements == offset) {
|
||||||
break
|
break
|
||||||
}
|
}
|
||||||
|
|
||||||
append(",\n")
|
append(",\n")
|
||||||
}
|
}
|
||||||
append("\n)")
|
append("\n)")
|
||||||
|
@ -182,7 +182,7 @@ internal class TestDoubleLinearOpsTensorAlgebra {
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
private inline fun DoubleLinearOpsTensorAlgebra.testSVDFor(tensor: DoubleTensor, epsilon: Double = 1e-10): Unit {
|
private fun DoubleLinearOpsTensorAlgebra.testSVDFor(tensor: DoubleTensor, epsilon: Double = 1e-10): Unit {
|
||||||
val svd = tensor.svd()
|
val svd = tensor.svd()
|
||||||
|
|
||||||
val tensorSVD = svd.first
|
val tensorSVD = svd.first
|
||||||
|
Loading…
Reference in New Issue
Block a user