102 lines
2.6 KiB
Kotlin
102 lines
2.6 KiB
Kotlin
// version 1.1.4-3
|
|
|
|
typealias Vector = DoubleArray
|
|
typealias Matrix = Array<Vector>
|
|
|
|
operator fun Matrix.times(other: Matrix): Matrix {
|
|
val rows1 = this.size
|
|
val cols1 = this[0].size
|
|
val rows2 = other.size
|
|
val cols2 = other[0].size
|
|
require(cols1 == rows2)
|
|
val result = Matrix(rows1) { Vector(cols2) }
|
|
for (i in 0 until rows1) {
|
|
for (j in 0 until cols2) {
|
|
for (k in 0 until rows2) {
|
|
result[i][j] += this[i][k] * other[k][j]
|
|
}
|
|
}
|
|
}
|
|
return result
|
|
}
|
|
|
|
fun pivotize(m: Matrix): Matrix {
|
|
val n = m.size
|
|
val im = Array(n) { Vector(n) }
|
|
for (i in 0 until n) im[i][i] = 1.0
|
|
for (i in 0 until n) {
|
|
var max = m[i][i]
|
|
var row = i
|
|
for (j in i until n) {
|
|
if (m[j][i] > max) {
|
|
max = m[j][i]
|
|
row = j
|
|
}
|
|
}
|
|
if (i != row) {
|
|
val t = im[i]
|
|
im[i] = im[row]
|
|
im[row] = t
|
|
}
|
|
}
|
|
return im
|
|
}
|
|
|
|
fun lu(a: Matrix): Array<Matrix> {
|
|
val n = a.size
|
|
val l = Array(n) { Vector(n) }
|
|
val u = Array(n) { Vector(n) }
|
|
val p = pivotize(a)
|
|
val a2 = p * a
|
|
|
|
for (j in 0 until n) {
|
|
l[j][j] = 1.0
|
|
for (i in 0 until j + 1) {
|
|
var sum = 0.0
|
|
for (k in 0 until i) sum += u[k][j] * l[i][k]
|
|
u[i][j] = a2[i][j] - sum
|
|
}
|
|
for (i in j until n) {
|
|
var sum2 = 0.0
|
|
for(k in 0 until j) sum2 += u[k][j] * l[i][k]
|
|
l[i][j] = (a2[i][j] - sum2) / u[j][j]
|
|
}
|
|
}
|
|
return arrayOf(l, u, p)
|
|
}
|
|
|
|
fun printMatrix(title: String, m: Matrix, f: String) {
|
|
val n = m.size
|
|
println("\n$title\n")
|
|
for (i in 0 until n) {
|
|
for (j in 0 until n) print("${f.format(m[i][j])} ")
|
|
println()
|
|
}
|
|
}
|
|
|
|
fun main(args: Array<String>) {
|
|
val a1 = arrayOf(
|
|
doubleArrayOf( 1.0, 3.0, 5.0),
|
|
doubleArrayOf( 2.0, 4.0, 7.0),
|
|
doubleArrayOf( 1.0, 1.0, 0.0)
|
|
)
|
|
val (l1, u1, p1) = lu(a1)
|
|
println("EXAMPLE 1:-")
|
|
printMatrix("A:", a1, "%1.0f")
|
|
printMatrix("L:", l1, "% 7.5f")
|
|
printMatrix("U:", u1, "% 8.5f")
|
|
printMatrix("P:", p1, "%1.0f")
|
|
|
|
val a2 = arrayOf(
|
|
doubleArrayOf(11.0, 9.0, 24.0, 2.0),
|
|
doubleArrayOf( 1.0, 5.0, 2.0, 6.0),
|
|
doubleArrayOf( 3.0, 17.0, 18.0, 1.0),
|
|
doubleArrayOf( 2.0, 5.0, 7.0, 1.0)
|
|
)
|
|
val (l2, u2, p2) = lu(a2)
|
|
println("\nEXAMPLE 2:-")
|
|
printMatrix("A:", a2, "%2.0f")
|
|
printMatrix("L:", l2, "%7.5f")
|
|
printMatrix("U:", u2, "%8.5f")
|
|
printMatrix("P:", p2, "%1.0f")
|
|
}
|