diff --git a/README.md b/README.md index 60bda66..cc61a3d 100644 --- a/README.md +++ b/README.md @@ -8,7 +8,7 @@ Another 15-459 Assignment(-ish). There are 2 versions with identical usage implemented with different libraries: - EJML version: No native dependencies, super lightweight, good for small input -- JBLAS version: Accelerated with BLAS, good for big input (see [jblas_dev](https://github.com/DEDZTBH/KuantumCircuitSim/tree/jblas_dev) branch) +- JBLAS version: Accelerated with BLAS, good for big input - JBLAS comes with default BLAS library for major OS's so it works right out of the box. For extreme performance, you can build JBLAS with custom BLAS library (like OpenBlas, ATLAS, cuBLAS) and put the library file(s) under library load path (for example, current working directory). See [JBLAS github page](https://github.com/jblas-project/jblas) for more detail. ### Table of Contents @@ -57,7 +57,7 @@ Options: --no_t, -q [false] -> Do not print circuit matrix in commandline after simulation if present --sequential, -s [false] -> Use sequential instead of concurrent implementation if present --init_state, -i [] -> Read custom initial joint state from csv if specified { String } - --binary_matrix, -b [false] -> Use binary format instead of csv for read/save circuit matrix if present + --binary_matrix, -b [false] -> Use binary format instead of csv for read/save circuit matrix if present (EJML version only) --help, -h -> Usage info ``` @@ -136,4 +136,4 @@ is stored in CSV file like 1,0,0,1 0,-1,1,0 ``` -Alternatively, circuit matrix can be saved/loaded in Java binary format (with -b option). It might be faster to save/load but it is not human readable. In case you want to read it, it is a `org.ejml.data.ZMatrixRMaj` (`org.jblas.ComplexDoubleMatrix` for JBLAS version). +(EJML version only) Alternatively, circuit matrix can be saved/loaded in Java binary format (with -b option). It might be faster to save/load but it is not human readable. In case you want to read it, it is a `org.ejml.data.ZMatrixRMaj`. diff --git a/build.gradle.kts b/build.gradle.kts index 9cfef3a..ba8215e 100644 --- a/build.gradle.kts +++ b/build.gradle.kts @@ -6,9 +6,13 @@ plugins { id("com.github.johnrengelman.shadow") version "6.0.0" } group = "com.dedztbh" -version = "1.4.3-ejml" +version = "1.4.4" + +val projectRoot = "${group}.kuantum" +val projectRootExclude = "/${projectRoot.replace('.', '/')}" val ejmlVersion = "0.39" +val jblasVersion = "1.2.5" repositories { mavenCentral() @@ -27,6 +31,7 @@ kotlin { dependencies { implementation("org.ejml:ejml-core:${ejmlVersion}") implementation("org.ejml:ejml-zdense:${ejmlVersion}") + implementation("org.jblas:jblas:${jblasVersion}") implementation("org.jetbrains.kotlinx:kotlinx-coroutines-core:1.3.9") implementation("com.github.cvb941:kotlin-parallel-operations:1.3") implementation("com.github.doyaaaaaken:kotlin-csv-jvm:0.11.1") @@ -44,18 +49,40 @@ tasks.withType { } tasks { - named("shadowJar") { - archiveBaseName.set("Kuantum") + register("shadowJarEjml", ShadowJar::class) { + archiveBaseName.set("Kuantum-ejml") + from(sourceSets.main.get().output) + configurations.add(project.configurations.runtime.get()) + configurations.add(project.configurations.implementation.get().apply { + isCanBeResolved = true + }) mergeServiceFiles() manifest { - attributes(mapOf("Main-Class" to "MainKt")) + attributes(mapOf("Main-Class" to "${projectRoot}.ejml.MainKt")) + } + dependencies { + exclude(dependency("org.jblas:jblas:${jblasVersion}")) } + exclude("${projectRootExclude}/jblas") minimize() } -} -tasks { - build { - dependsOn(shadowJar) + register("shadowJarJBLAS", ShadowJar::class) { + archiveBaseName.set("Kuantum-jblas") + from(sourceSets.main.get().output) + configurations.add(project.configurations.runtime.get()) + configurations.add(project.configurations.implementation.get().apply { + isCanBeResolved = true + }) + mergeServiceFiles() + manifest { + attributes(mapOf("Main-Class" to "${projectRoot}.jblas.MainKt")) + } + dependencies { + exclude(dependency("org.ejml:ejml-core:${ejmlVersion}")) + exclude(dependency("org.ejml:ejml-zdense:${ejmlVersion}")) + } + exclude("${projectRootExclude}/ejml") + minimize() } } \ No newline at end of file diff --git a/src/main/kotlin/Config.kt b/src/main/kotlin/com/dedztbh/kuantum/common/Config.kt similarity index 92% rename from src/main/kotlin/Config.kt rename to src/main/kotlin/com/dedztbh/kuantum/common/Config.kt index a0a37d6..e014b9a 100644 --- a/src/main/kotlin/Config.kt +++ b/src/main/kotlin/com/dedztbh/kuantum/common/Config.kt @@ -1,9 +1,10 @@ +package com.dedztbh.kuantum.common + import kotlinx.cli.ArgParser import kotlinx.cli.ArgType import kotlinx.cli.default import kotlinx.cli.optional - /** * Created by DEDZTBH on 2020/09/27. * Project KuantumCircuitSim @@ -11,6 +12,8 @@ import kotlinx.cli.optional const val DEFAULT_N = 5 +const val CMD = "java -jar Kuantum.jar" + class Config(parser: ArgParser) { val input by parser.argument(ArgType.String, description = "Input file") val operator by parser.argument(ArgType.String, description = "Operator name") @@ -44,6 +47,6 @@ class Config(parser: ArgParser) { val binary_matrix by parser.option( ArgType.Boolean, shortName = "b", - description = "Use binary format instead of csv for read/save circuit matrix if present" + description = "Use binary format instead of csv for read/save circuit matrix if present (EJML version only)" ).default(false) } \ No newline at end of file diff --git a/src/main/kotlin/operator/Operator.kt b/src/main/kotlin/com/dedztbh/kuantum/common/Operator.kt similarity index 57% rename from src/main/kotlin/operator/Operator.kt rename to src/main/kotlin/com/dedztbh/kuantum/common/Operator.kt index 8735c46..69e7d47 100644 --- a/src/main/kotlin/operator/Operator.kt +++ b/src/main/kotlin/com/dedztbh/kuantum/common/Operator.kt @@ -1,6 +1,5 @@ -package operator +package com.dedztbh.kuantum.common -import Config import kotlinx.coroutines.CoroutineScope /** @@ -9,7 +8,7 @@ import kotlinx.coroutines.CoroutineScope */ /** - * Children class must have a constructor that takes an Config object and a coroutine scope + * Children class must have a constructor that takes an common.Config object and a coroutine scope */ interface Operator { suspend fun runCmd(cmd: String): Int @@ -18,8 +17,8 @@ interface Operator { companion object { @JvmStatic - fun get(config: Config, scope: CoroutineScope) = - Class.forName("operator.${if (config.sequential) "" else "P"}${config.operator}") + fun get(config: Config, scope: CoroutineScope, lib: String) = + Class.forName("com.dedztbh.kuantum.${lib}.operator.${if (config.sequential) "" else "P"}${config.operator}") .getDeclaredConstructor(Config::class.java, CoroutineScope::class.java) .newInstance(config, scope) as Operator } diff --git a/src/main/kotlin/main.kt b/src/main/kotlin/com/dedztbh/kuantum/common/main.kt similarity index 77% rename from src/main/kotlin/main.kt rename to src/main/kotlin/com/dedztbh/kuantum/common/main.kt index 83ffae9..f41f722 100644 --- a/src/main/kotlin/main.kt +++ b/src/main/kotlin/com/dedztbh/kuantum/common/main.kt @@ -1,7 +1,8 @@ +package com.dedztbh.kuantum.common + import kotlinx.cli.ArgParser import kotlinx.coroutines.Dispatchers import kotlinx.coroutines.runBlocking -import operator.Operator import java.io.File /** @@ -9,15 +10,15 @@ import java.io.File * Project KuantumCircuitSim */ -fun main(args: Array) { - val parser = ArgParser("java -jar Kuantum.jar") +fun main(args: Array, lib: String) { + val parser = ArgParser(CMD) val config = Config(parser) parser.parse(args) reader = File(config.input).bufferedReader() runBlocking(Dispatchers.Default) { - val operator = Operator.get(config, this) + val operator = Operator.get(config, this, lib) var cmd = read() while (true) { diff --git a/src/main/kotlin/util.kt b/src/main/kotlin/com/dedztbh/kuantum/common/util.kt similarity index 96% rename from src/main/kotlin/util.kt rename to src/main/kotlin/com/dedztbh/kuantum/common/util.kt index 620b9bc..8f115e6 100644 --- a/src/main/kotlin/util.kt +++ b/src/main/kotlin/com/dedztbh/kuantum/common/util.kt @@ -1,3 +1,5 @@ +package com.dedztbh.kuantum.common + import java.io.BufferedReader import java.util.* import kotlin.math.pow diff --git a/src/main/kotlin/com/dedztbh/kuantum/ejml/main.kt b/src/main/kotlin/com/dedztbh/kuantum/ejml/main.kt new file mode 100644 index 0000000..ae67312 --- /dev/null +++ b/src/main/kotlin/com/dedztbh/kuantum/ejml/main.kt @@ -0,0 +1,10 @@ +package com.dedztbh.kuantum.ejml + +import com.dedztbh.kuantum.common.main + +/** + * Created by DEDZTBH on 2020/10/11. + * Project KuantumCircuitSim + */ + +fun main(args: Array) = main(args, "ejml") \ No newline at end of file diff --git a/src/main/kotlin/com/dedztbh/kuantum/ejml/matrix/CMatrixIO.kt b/src/main/kotlin/com/dedztbh/kuantum/ejml/matrix/CMatrixIO.kt new file mode 100644 index 0000000..d03fb61 --- /dev/null +++ b/src/main/kotlin/com/dedztbh/kuantum/ejml/matrix/CMatrixIO.kt @@ -0,0 +1,112 @@ +package com.dedztbh.kuantum.ejml.matrix + +import com.github.doyaaaaaken.kotlincsv.client.CsvReader +import com.github.doyaaaaaken.kotlincsv.client.CsvWriter +import org.ejml.UtilEjml +import org.ejml.data.MatrixSparse +import org.ejml.data.MatrixType +import org.ejml.ops.MatrixIO +import java.io.* +import java.text.DecimalFormat + +/** + * Created by DEDZTBH on 2020/09/27. + * Project KuantumCircuitSim + */ + +fun getMatrixType(mat: CMatrix): String { + return if (mat.type == MatrixType.UNSPECIFIED) { + mat.javaClass.simpleName + } else { + mat.type.name + } +} + +fun printTypeSize(out: PrintStream, mat: CMatrix) { + if (mat is MatrixSparse) { + val m = mat as MatrixSparse + out.println( + "Type = " + getMatrixType(mat) + " , rows = " + mat.getNumRows() + + " , cols = " + mat.getNumCols() + " , nz_length = " + m.nonZeroLength + ) + } else { + out.println("Type = " + getMatrixType(mat) + " , rows = " + mat.getNumRows() + " , cols = " + mat.getNumCols()) + } +} + +fun padSpace(builder: java.lang.StringBuilder, length: Int): String { + builder.delete(0, builder.length) + for (i in 0 until length) { + builder.append(' ') + } + return builder.toString() +} + +fun saveCsv(A: CMatrix, fileName: String) = + CsvWriter().writeAll(A.data.asList().windowed(A.rowStride, A.rowStride), fileName) + + +fun loadCsv(fileName: String): CMatrix = + CMatrix(CsvReader().open(fileName) { + readAllAsSequence().map { + DoubleArray(it.size) { i -> it[i].toDouble() } + }.toMutableList().toTypedArray() + }) + +fun saveBin(A: T, fileName: String) { + val fileStream = FileOutputStream(fileName) + val stream = ObjectOutputStream(fileStream) + try { + stream.writeObject(A) + stream.flush() + } finally { + // clean up + fileStream.use { + stream.close() + } + } +} + +@Suppress("UNCHECKED_CAST") +fun loadBin(fileName: String): T { + val fileStream = FileInputStream(fileName) + val stream = ObjectInputStream(fileStream) + val ret: T + try { + ret = stream.readObject() as T + if (stream.available() != 0) { + throw RuntimeException("File not completely read?") + } + } catch (e: ClassNotFoundException) { + throw RuntimeException(e) + } + stream.close() + return ret +} + +fun CMatrix.printFancy2( + allssket: List, + out: PrintStream = System.out, + length: Int = MatrixIO.DEFAULT_LENGTH, +) { +// printTypeSize(out, this) + val format = DecimalFormat("#") + val builder = StringBuilder(length) + val cols = numCols + val c = CNum() + var i = 0 + for (y in 0 until numRows) { + for (x in 0 until cols) { + get(y, x, c) + var real = UtilEjml.fancyString(c.real, format, length, 4) + var img = UtilEjml.fancyString(c.imaginary, format, length, 4) + real += padSpace(builder, length - real.length) + img = img + "i" + padSpace(builder, length - img.length) + out.print("${allssket[i++]}: $real + $img") + if (x < numCols - 1) { + out.print(" , ") + } + } + out.println() + } +} \ No newline at end of file diff --git a/src/main/kotlin/matrix/const.kt b/src/main/kotlin/com/dedztbh/kuantum/ejml/matrix/const.kt similarity index 97% rename from src/main/kotlin/matrix/const.kt rename to src/main/kotlin/com/dedztbh/kuantum/ejml/matrix/const.kt index a4bd7f1..d3b43c1 100644 --- a/src/main/kotlin/matrix/const.kt +++ b/src/main/kotlin/com/dedztbh/kuantum/ejml/matrix/const.kt @@ -1,4 +1,4 @@ -package matrix +package com.dedztbh.kuantum.ejml.matrix /** * Created by DEDZTBH on 2020/09/25. diff --git a/src/main/kotlin/matrix/op.kt b/src/main/kotlin/com/dedztbh/kuantum/ejml/matrix/op.kt similarity index 96% rename from src/main/kotlin/matrix/op.kt rename to src/main/kotlin/com/dedztbh/kuantum/ejml/matrix/op.kt index 0ca90b0..b2d3822 100644 --- a/src/main/kotlin/matrix/op.kt +++ b/src/main/kotlin/com/dedztbh/kuantum/ejml/matrix/op.kt @@ -1,5 +1,4 @@ -package matrix - +package com.dedztbh.kuantum.ejml.matrix /** * Created by DEDZTBH on 2020/09/22. diff --git a/src/main/kotlin/matrix/type.kt b/src/main/kotlin/com/dedztbh/kuantum/ejml/matrix/type.kt similarity index 89% rename from src/main/kotlin/matrix/type.kt rename to src/main/kotlin/com/dedztbh/kuantum/ejml/matrix/type.kt index 21a1b15..c0a20e9 100644 --- a/src/main/kotlin/matrix/type.kt +++ b/src/main/kotlin/com/dedztbh/kuantum/ejml/matrix/type.kt @@ -1,4 +1,4 @@ -package matrix +package com.dedztbh.kuantum.ejml.matrix import org.ejml.data.Complex_F64 import org.ejml.data.ZMatrixRMaj diff --git a/src/main/kotlin/operator/AllInit.kt b/src/main/kotlin/com/dedztbh/kuantum/ejml/operator/AllInit.kt similarity index 67% rename from src/main/kotlin/operator/AllInit.kt rename to src/main/kotlin/com/dedztbh/kuantum/ejml/operator/AllInit.kt index 924860e..82cfb69 100644 --- a/src/main/kotlin/operator/AllInit.kt +++ b/src/main/kotlin/com/dedztbh/kuantum/ejml/operator/AllInit.kt @@ -1,12 +1,10 @@ -package operator +package com.dedztbh.kuantum.ejml.operator -import Config +import com.dedztbh.kuantum.common.Config +import com.dedztbh.kuantum.ejml.matrix.* import com.lukaskusik.coroutines.transformations.map.mapParallel -import com.lukaskusik.coroutines.transformations.reduce.reduceParallel import kotlinx.coroutines.CoroutineScope import kotlinx.coroutines.async -import matrix.* -import matrix.CMatrixIO.printFancy2 /** * Created by DEDZTBH on 2020/09/22. @@ -19,8 +17,9 @@ class AllInit(config: Config, scope: CoroutineScope) : TFinder(config, scope) { println("\nFinal states: ") alls.forEachIndexed { idx, arr -> println("Init ${allssket[idx]}") - var jointState = I1 - arr.forEach { i -> + val it = arr.iterator() + var jointState = if (it.next() == 0) KET0 else KET1 + it.forEachRemaining { i -> jointState = jointState kron (if (i == 0) KET0 else KET1) } (opMatrix * jointState).printFancy2(allssket = allssket) @@ -35,9 +34,12 @@ class PAllInit(config: Config, scope: CoroutineScope) : PTFinder(config, scope) println("\nFinal states: ") alls.mapParallel { arr -> scope.async { - opMatrix * - arr.mapParallel { if (it == 0) KET0 else KET1 } - .reduceParallel { a, b -> a kron b } + val it = arr.iterator() + var jointState = if (it.next() == 0) KET0 else KET1 + it.forEachRemaining { i -> + jointState = jointState kron (if (i == 0) KET0 else KET1) + } + opMatrix * jointState } }.forEachIndexed { i, mat -> println("Init ${allssket[i]}") diff --git a/src/main/kotlin/operator/TFinder.kt b/src/main/kotlin/com/dedztbh/kuantum/ejml/operator/TFinder.kt similarity index 96% rename from src/main/kotlin/operator/TFinder.kt rename to src/main/kotlin/com/dedztbh/kuantum/ejml/operator/TFinder.kt index 84b6282..e0b1f9b 100644 --- a/src/main/kotlin/operator/TFinder.kt +++ b/src/main/kotlin/com/dedztbh/kuantum/ejml/operator/TFinder.kt @@ -1,13 +1,10 @@ -package operator +package com.dedztbh.kuantum.ejml.operator -import Config -import allStates +import com.dedztbh.kuantum.common.* +import com.dedztbh.kuantum.ejml.matrix.* import com.lukaskusik.coroutines.transformations.reduce.reduceParallel import kotlinx.coroutines.CoroutineScope import kotlinx.coroutines.async -import matrix.* -import readDouble -import readInt import java.util.concurrent.ConcurrentHashMap import kotlin.math.PI import kotlin.math.cos @@ -42,9 +39,9 @@ open class TFinder(val config: Config, val scope: CoroutineScope) : Operator { @JvmField var opMatrix = if (config.input_matrix.isNotBlank()) { if (config.binary_matrix) - CMatrixIO.loadBin(config.input_matrix) + loadBin(config.input_matrix) else - CMatrixIO.loadCsv(config.input_matrix) + loadCsv(config.input_matrix) } else IN2 @JvmField @@ -186,9 +183,9 @@ open class TFinder(val config: Config, val scope: CoroutineScope) : Operator { override suspend fun done() { if (config.output.isNotBlank()) { if (config.binary_matrix) - CMatrixIO.saveBin(opMatrix, config.output) + saveBin(opMatrix, config.output) else - CMatrixIO.saveCsv(opMatrix, config.output) + saveCsv(opMatrix, config.output) } } diff --git a/src/main/kotlin/operator/Tester.kt b/src/main/kotlin/com/dedztbh/kuantum/ejml/operator/Tester.kt similarity index 96% rename from src/main/kotlin/operator/Tester.kt rename to src/main/kotlin/com/dedztbh/kuantum/ejml/operator/Tester.kt index e21d034..1381bce 100644 --- a/src/main/kotlin/operator/Tester.kt +++ b/src/main/kotlin/com/dedztbh/kuantum/ejml/operator/Tester.kt @@ -1,13 +1,12 @@ -package operator +package com.dedztbh.kuantum.ejml.operator -import Config +import com.dedztbh.kuantum.common.Config +import com.dedztbh.kuantum.common.readInt +import com.dedztbh.kuantum.common.upperBound +import com.dedztbh.kuantum.ejml.matrix.* import kotlinx.coroutines.CoroutineScope import kotlinx.coroutines.GlobalScope import kotlinx.coroutines.async -import matrix.* -import matrix.CMatrixIO.printFancy2 -import readInt -import upperBound import kotlin.math.sqrt import kotlin.random.Random @@ -21,7 +20,7 @@ fun getJointState(initState: String, jointStateSize: Int) = CMatrix(jointStateSize, 1).apply { set(0, 0, 1.0, 0.0) } - else CMatrixIO.loadCsv(initState) + else loadCsv(initState) class Tester(config: Config, scope: CoroutineScope) : TFinder(config, scope) { /** 2^N by 1 column vector */ diff --git a/src/main/kotlin/com/dedztbh/kuantum/jblas/main.kt b/src/main/kotlin/com/dedztbh/kuantum/jblas/main.kt new file mode 100644 index 0000000..35e6b9a --- /dev/null +++ b/src/main/kotlin/com/dedztbh/kuantum/jblas/main.kt @@ -0,0 +1,14 @@ +package com.dedztbh.kuantum.jblas + +import com.dedztbh.kuantum.common.main +import org.jblas.util.Logger + +/** + * Created by DEDZTBH on 2020/10/11. + * Project KuantumCircuitSim + */ + +fun main(args: Array) { + Logger.getLogger().setLevel(Logger.WARNING) + main(args, "jblas") +} \ No newline at end of file diff --git a/src/main/kotlin/com/dedztbh/kuantum/jblas/matrix/CMatrixIO.kt b/src/main/kotlin/com/dedztbh/kuantum/jblas/matrix/CMatrixIO.kt new file mode 100644 index 0000000..b916fba --- /dev/null +++ b/src/main/kotlin/com/dedztbh/kuantum/jblas/matrix/CMatrixIO.kt @@ -0,0 +1,53 @@ +package com.dedztbh.kuantum.jblas.matrix + +import com.github.doyaaaaaken.kotlincsv.client.CsvReader +import com.github.doyaaaaaken.kotlincsv.client.CsvWriter +import java.io.PrintStream +import java.io.PrintWriter + +/** + * Created by DEDZTBH on 2020/09/27. + * Project KuantumCircuitSim + */ + +fun saveCsv(A: CMatrix, fileName: String) { + val rowStride = A.columns * 2 + CsvWriter().writeAll(A.transpose().data.asList().windowed(rowStride, rowStride), fileName) +} + +fun loadCsv(fileName: String): CMatrix = + CsvReader().open(fileName) { + val m = readAllAsSequence().toList() + CMatrix(m[0].size / 2, m.size).apply { + data = m.flatMap { it.map { it.toDouble() } }.toDoubleArray() + }.transpose() + } + +fun CMatrix.toStringFancy( + allssket: List? = null +): String { + val s = StringBuilder() + for (i in 0 until rows) { + allssket?.let { s.append("${it[i]}: ") } + for (j in 0 until columns) { + s.append(get(i, j)) + if (j < columns - 1) s.append(", ") + } + if (i < rows - 1) s.append("\n") + } + return s.toString() +} + +fun CMatrix.printFancy( + allssket: List? = null, + o: PrintStream = System.out +) = PrintWriter(o, false).apply { + for (i in 0 until rows) { + allssket?.let { print("${it[i]}: ") } + for (j in 0 until columns) { + print(get(i, j)) + if (j < columns - 1) print(", ") + } + if (i < rows - 1) print("\n") + } +}.flush() \ No newline at end of file diff --git a/src/main/kotlin/com/dedztbh/kuantum/jblas/matrix/const.kt b/src/main/kotlin/com/dedztbh/kuantum/jblas/matrix/const.kt new file mode 100644 index 0000000..f172122 --- /dev/null +++ b/src/main/kotlin/com/dedztbh/kuantum/jblas/matrix/const.kt @@ -0,0 +1,67 @@ +package com.dedztbh.kuantum.jblas.matrix + +import org.jblas.ComplexDoubleMatrix +import org.jblas.DoubleMatrix + +/** + * Created by DEDZTBH on 2020/09/25. + * Project KuantumCircuitSim + */ + +/** 2^(-1/2) */ +const val HALF_AMPL = 0.70710678118654757273731092936941422522068023681640625 + +/** Don't change these constant matrices! */ +val NOT: ComplexDoubleMatrix = CMatrix( + 2, 2, + 0.0, 0.0, 1.0, 0.0, + 1.0, 0.0, 0.0, 0.0 +).transpose() +val H: ComplexDoubleMatrix = CMatrix( + 2, 2, + HALF_AMPL, 0.0, HALF_AMPL, 0.0, + HALF_AMPL, 0.0, -HALF_AMPL, 0.0 +).transpose() +val I1 = CMatrix(DoubleMatrix.eye(1), null) +val I2 = CMatrix(DoubleMatrix.eye(2), null) + +val KET0: ComplexDoubleMatrix = CMatrix( + 2, 1, + 1.0, 0.0, + 0.0, 0.0 +).transpose() +val KET1: ComplexDoubleMatrix = CMatrix( + 2, 1, + 0.0, 0.0, + 1.0, 0.0 +).transpose() + +val KETBRA0: ComplexDoubleMatrix = CMatrix.diag(CMatrix(doubleArrayOf(1.0, 0.0, 0.0, 0.0)).transpose()) +val KETBRA1: ComplexDoubleMatrix = CMatrix.diag(CMatrix(doubleArrayOf(0.0, 0.0, 1.0, 0.0)).transpose()) +val SQRT_NOT: ComplexDoubleMatrix = CMatrix( + 2, 2, + 0.5, 0.5, 0.5, -0.5, + 0.5, -0.5, 0.5, 0.5 +).transpose() +val SQRT_NOT_DAG: ComplexDoubleMatrix = SQRT_NOT.hermitian() +val Y: ComplexDoubleMatrix = CMatrix( + 2, 2, + 0.0, 0.0, 0.0, -1.0, + 0.0, 1.0, 0.0, 0.0 +).transpose() +val Z: ComplexDoubleMatrix = CMatrix( + 2, 2, + 1.0, 0.0, 0.0, 0.0, + 0.0, 0.0, -1.0, 0.0 +).transpose() +val S: ComplexDoubleMatrix = CMatrix( + 2, 2, + 1.0, 0.0, 0.0, 0.0, + 0.0, 0.0, 0.0, 1.0 +).transpose() +val T: ComplexDoubleMatrix = CMatrix( + 2, 2, + 1.0, 0.0, 0.0, 0.0, + 0.0, 0.0, HALF_AMPL, HALF_AMPL +).transpose() +val TDag: ComplexDoubleMatrix = T.hermitian() diff --git a/src/main/kotlin/com/dedztbh/kuantum/jblas/matrix/op.kt b/src/main/kotlin/com/dedztbh/kuantum/jblas/matrix/op.kt new file mode 100644 index 0000000..168e440 --- /dev/null +++ b/src/main/kotlin/com/dedztbh/kuantum/jblas/matrix/op.kt @@ -0,0 +1,53 @@ +package com.dedztbh.kuantum.jblas.matrix + +import org.jblas.ComplexDoubleMatrix + + +/** + * Created by DEDZTBH on 2020/09/22. + * Project KuantumCircuitSim + */ + +/** + * Kronecker product of two complex matrices + */ + +infix fun CMatrix.kron(B: CMatrix) = CMatrix(rows * B.rows, columns * B.columns).also { C -> + var jbcol = 0 + for (j in 0 until columns) { + // jbcol = j * B.columns + var ibrow = 0 + for (i in 0 until rows) { + // ibrow = i * B.rows + val a = get(i, j) + var c = jbcol + for (colB in 0 until B.columns) { + // c = jbcol + colB + var r = ibrow + for (rowB in 0 until B.rows) { + // r = ibrow + rowB + C.put(r++, c, a * (B.get(rowB, colB))) + } + ++c + } + ibrow += B.rows + } + jbcol += B.columns + } +} + +/** + * Convenient pure ops + */ +inline operator fun CMatrix.plus(B: CMatrix): ComplexDoubleMatrix = add(B) + +inline operator fun CMatrix.times(B: CMatrix): ComplexDoubleMatrix = mmul(B) + +val CNum.magnitude2 + get() = let { + val r = real() + val i = imag() + r * r + i * i + } + +inline operator fun CNum.times(B: CNum): CNum = mul(B) \ No newline at end of file diff --git a/src/main/kotlin/com/dedztbh/kuantum/jblas/matrix/type.kt b/src/main/kotlin/com/dedztbh/kuantum/jblas/matrix/type.kt new file mode 100644 index 0000000..29a1090 --- /dev/null +++ b/src/main/kotlin/com/dedztbh/kuantum/jblas/matrix/type.kt @@ -0,0 +1,12 @@ +package com.dedztbh.kuantum.jblas.matrix + +import org.jblas.ComplexDouble +import org.jblas.ComplexDoubleMatrix + +/** + * Created by DEDZTBH on 2020/09/25. + * Project KuantumCircuitSim + */ + +typealias CMatrix = ComplexDoubleMatrix +typealias CNum = ComplexDouble diff --git a/src/main/kotlin/com/dedztbh/kuantum/jblas/operator/AllInit.kt b/src/main/kotlin/com/dedztbh/kuantum/jblas/operator/AllInit.kt new file mode 100644 index 0000000..8235adf --- /dev/null +++ b/src/main/kotlin/com/dedztbh/kuantum/jblas/operator/AllInit.kt @@ -0,0 +1,50 @@ +package com.dedztbh.kuantum.jblas.operator + +import com.dedztbh.kuantum.common.Config +import com.dedztbh.kuantum.jblas.matrix.* +import com.lukaskusik.coroutines.transformations.map.mapParallel +import kotlinx.coroutines.CoroutineScope +import kotlinx.coroutines.async + +/** + * Created by DEDZTBH on 2020/09/22. + * Project KuantumCircuitSim + */ + +class AllInit(config: Config, scope: CoroutineScope) : TFinder(config, scope) { + override suspend fun printResult() { + super.printResult() + println("\nFinal states: ") + alls.forEachIndexed { idx, arr -> + println("Init ${allssket[idx]}") + val it = arr.iterator() + var jointState = if (it.next() == 0) KET0 else KET1 + it.forEachRemaining { i -> + jointState = jointState kron (if (i == 0) KET0 else KET1) + } + (opMatrix * jointState).printFancy(allssket) + println() + } + } +} + +class PAllInit(config: Config, scope: CoroutineScope) : PTFinder(config, scope) { + override suspend fun printResult() { + super.printResult() + println("\nFinal states: ") + alls.mapParallel { arr -> + scope.async { + val it = arr.iterator() + var jointState = if (it.next() == 0) KET0 else KET1 + it.forEachRemaining { i -> + jointState = jointState kron (if (i == 0) KET0 else KET1) + } + opMatrix * jointState + } + }.forEachIndexed { i, mat -> + println("Init ${allssket[i]}") + mat.await().printFancy(allssket) + println() + } + } +} \ No newline at end of file diff --git a/src/main/kotlin/com/dedztbh/kuantum/jblas/operator/TFinder.kt b/src/main/kotlin/com/dedztbh/kuantum/jblas/operator/TFinder.kt new file mode 100644 index 0000000..fdc3ea1 --- /dev/null +++ b/src/main/kotlin/com/dedztbh/kuantum/jblas/operator/TFinder.kt @@ -0,0 +1,280 @@ +package com.dedztbh.kuantum.jblas.operator + +import com.dedztbh.kuantum.common.* +import com.dedztbh.kuantum.jblas.matrix.* +import com.lukaskusik.coroutines.transformations.reduce.reduceParallel +import kotlinx.coroutines.CoroutineScope +import kotlinx.coroutines.async +import org.jblas.DoubleMatrix +import java.util.concurrent.ConcurrentHashMap +import kotlin.math.PI +import kotlin.math.cos +import kotlin.math.sin + +/** + * Created by DEDZTBH on 2020/09/22. + * Project KuantumCircuitSim + */ + +open class TFinder(val config: Config, val scope: CoroutineScope) : Operator { + @JvmField + val N = config.N + + /** all possible states (left-to-right) */ + @JvmField + val alls = allStates(N) + + /** all possible states represented with ket (right-to-left) */ + @JvmField + val allss = alls.map { arr -> arr.joinToString("") } + + @JvmField + val allssket = allss.map { "|$it>" } + + @JvmField + val jointStateSize = alls.size + + @JvmField + val IN2 = CMatrix(DoubleMatrix.eye(jointStateSize), null) + + @JvmField + var opMatrix = if (config.input_matrix.isNotBlank()) { + loadCsv(config.input_matrix) + } else IN2 + + @JvmField + val IKronTable = (0..N).scan(I1) { it, _ -> it kron I2 }.toTypedArray() + + open fun getHashMap(): MutableMap = HashMap() + + @JvmField + val matrix0CtrlCache: Array> = Array(N) { getHashMap() } + fun get0CtrlMatrix(i: Int, mat: CMatrix, cache: Boolean = true): CMatrix { + if (cache) matrix0CtrlCache[i][mat]?.let { return it } + val res = IKronTable[i] kron mat kron IKronTable[N - i - 1] + if (cache) matrix0CtrlCache[i][mat] = res + return res + } + + /** Control matrix generation based on this great article: + * http://www.sakkaris.com/tutorials/quantum_control_gates.html */ + @JvmField + val matrix1CtrlCache: Array>> = Array(N) { Array(N) { getHashMap() } } + fun get1CtrlMatrix(i: Int, j: Int, mat: CMatrix): CMatrix { + matrix1CtrlCache[i][j][mat]?.let { return it } + val res = get0CtrlMatrix(i, KETBRA0) + when { + i < j -> IKronTable[i] kron KETBRA1 kron + IKronTable[j - i - 1] kron + mat kron IKronTable[N - j - 1] + i > j -> IKronTable[j] kron mat kron + IKronTable[i - j - 1] kron + KETBRA1 kron IKronTable[N - i - 1] + else -> throw IllegalArgumentException("Control qubit is same as affected qubit") + } + matrix1CtrlCache[i][j][mat] = res + return res + } + + @JvmField + val ccnotCache: Array>> = Array(N) { Array(N) { getHashMap() } } + fun getCCNotMatrix(i: Int, j: Int, k: Int): CMatrix { + /** using Sleator-Weinfurter construction */ + ccnotCache[i][j][k]?.let { return it } + val cnotij = get1CtrlMatrix(i, j, NOT) + val res = get1CtrlMatrix(i, k, SQRT_NOT) * + cnotij * + get1CtrlMatrix(j, k, SQRT_NOT_DAG) * + cnotij * + get1CtrlMatrix(j, k, SQRT_NOT) + ccnotCache[i][j][k] = res + return res + /** Alternative implementation: + * https://en.wikipedia.org/wiki/Toffoli_gate + val Hk = get0CtrlMatrix(k, H) + val TDagj = get0CtrlMatrix(j, TDag) + val TDagk = get0CtrlMatrix(k, TDag) + val Ti = get0CtrlMatrix(i, T) + val Tj = get0CtrlMatrix(j, T) + val Tk = get0CtrlMatrix(k, T) + val CNotij = get1CtrlMatrix(i, j, NOT) + val CNotik = get1CtrlMatrix(i, k, NOT) + val CNotjk = get1CtrlMatrix(j, k, NOT) + return CNotij * TDagj * Ti * + Hk * CNotij * Tk * Tj * + CNotik * TDagk * CNotjk * + Tk * CNotik * TDagk * CNotjk * Hk */ + } + + @JvmField + val rotCache: Array> = Array(N) { getHashMap() } + fun getRotMatrix(i: Int, d: Double): CMatrix { + rotCache[i][d]?.let { return it } + val rad = d * PI / 180.0 + val sine = sin(rad) + val cosine = cos(rad) + + /** This is column major! + * It's actually + * cos -sin + * sin cos */ + val rotMat = CMatrix( + 2, 2, + cosine, 0.0, sine, 0.0, + -sine, 0.0, cosine, 0.0 + ) + return get0CtrlMatrix(i, rotMat, false).also { + rotCache[i][d] = it + } + } + + @JvmField + val swapCache: Array> = Array(N) { getHashMap() } + fun getSwapMatrix(i: Int, j: Int): CMatrix { + swapCache[i][j]?.let { return it } + /** https://algassert.com/post/1717 + * Swap implemented with 3 CNots */ + val cnot0 = get1CtrlMatrix(i, j, NOT) + return (cnot0 * get1CtrlMatrix(j, i, NOT) * cnot0).also { + swapCache[i][j] = it + } + } + + @JvmField + val cswapCache: Array>> = Array(N) { Array(N) { getHashMap() } } + fun getCSwapMatrix(i: Int, j: Int, k: Int): CMatrix { + cswapCache[i][j][k]?.let { return it } + /** https://quantumcomputing.stackexchange.com/ + * questions/9342/how-to-implement-a-fredkin + * -gate-using-toffoli-and-cnots */ + val cnotkj = get1CtrlMatrix(k, j, NOT) + return (cnotkj * getCCNotMatrix(i, j, k) * cnotkj).also { + cswapCache[i][j][k] = it + } + } + + override suspend fun runCmd(cmd: String): Int { + val i = readInt() + val newOp = when (cmd) { + "NOT" -> get0CtrlMatrix(i, NOT) + "HADAMARD", "H" -> get0CtrlMatrix(i, H) + "CNOT" -> get1CtrlMatrix(i, readInt(), NOT) + "SWAP" -> getSwapMatrix(i, readInt()) + "CCNOT" -> getCCNotMatrix(i, readInt(), readInt()) + "CSWAP" -> getCSwapMatrix(i, readInt(), readInt()) + "Y" -> get0CtrlMatrix(i, Y) + "Z" -> get0CtrlMatrix(i, Z) + "S" -> get0CtrlMatrix(i, S) + "T" -> get0CtrlMatrix(i, T) + "TDAG" -> get0CtrlMatrix(i, TDag) + "SQRTNOT" -> get0CtrlMatrix(i, SQRT_NOT) + "SQRTNOTDAG" -> get0CtrlMatrix(i, SQRT_NOT_DAG) +// "SQRTSWAP" -> { +// TODO: Implement SqrtSwap +// } + "ROT" -> getRotMatrix(i, readDouble()) + "CZ" -> get1CtrlMatrix(i, readInt(), Z) + else -> { + System.err.println("Unknown command \"${cmd}\". Stop reading commands.") + return -1 + } + } + opMatrix = newOp * opMatrix + return 0 + } + + override suspend fun done() { + if (config.output.isNotBlank()) { + saveCsv(opMatrix, config.output) + } + } + + override suspend fun printResult() { + if (!config.no_t) { + println("Transformation: ") + opMatrix.printFancy() + println() + } + } +} + +const val CONCURRENT_MATRIX = 4096 + +open class PTFinder(config: Config, scope: CoroutineScope) : TFinder(config, scope) { + + @JvmField + var reversedNewOps = mutableListOf(scope.async { opMatrix }) + + override fun getHashMap(): MutableMap = ConcurrentHashMap() + + override suspend fun runCmd(cmd: String): Int { + val i = readInt() + scope.run { + reversedNewOps.add(when (cmd) { + "NOT" -> async { get0CtrlMatrix(i, NOT) } + "HADAMARD", "H" -> async { get0CtrlMatrix(i, H) } + "CNOT" -> { + val j = readInt() + async { get1CtrlMatrix(i, j, NOT) } + } + "SWAP" -> { + /** https://algassert.com/post/1717 + * Swap implemented with 3 CNots */ + val j = readInt() + async { getSwapMatrix(i, j) } + } + "CCNOT" -> { + val j = readInt() + val k = readInt() + async { getCCNotMatrix(i, j, k) } + } + "CSWAP" -> { + /** https://quantumcomputing.stackexchange.com/ + * questions/9342/how-to-implement-a-fredkin + * -gate-using-toffoli-and-cnots */ + val j = readInt() + val k = readInt() + async { getCSwapMatrix(i, j, k) } + } + "Y" -> async { get0CtrlMatrix(i, Y) } + "Z" -> async { get0CtrlMatrix(i, Z) } + "S" -> async { get0CtrlMatrix(i, S) } + "T" -> async { get0CtrlMatrix(i, T) } + "TDAG" -> async { get0CtrlMatrix(i, TDag) } + "SQRTNOT" -> async { get0CtrlMatrix(i, SQRT_NOT) } + "SQRTNOTDAG" -> async { get0CtrlMatrix(i, SQRT_NOT_DAG) } +// "SQRTSWAP" -> { +// TODO: Implement SqrtSwap +// } + "ROT" -> { + val deg = readDouble() + async { getRotMatrix(i, deg) } + } + "CZ" -> { + val j = readInt() + async { get1CtrlMatrix(i, j, Z) } + } + else -> { + System.err.println("Unknown command \"${cmd}\". Stop reading commands.") + return -1 + } + }) + if (reversedNewOps.size >= CONCURRENT_MATRIX) { + // We have enough matrices, do a reduction to prevent OOM + val newNewOps = reduceOps() + reversedNewOps = mutableListOf(async { newNewOps }) + } + } + return 0 + } + + suspend fun reduceOps() = + reversedNewOps.map { it.await() } + .reduceParallel { d1, d2 -> + d2 * d1 + } + + override suspend fun done() { + opMatrix = reduceOps() + super.done() + } +} \ No newline at end of file diff --git a/src/main/kotlin/com/dedztbh/kuantum/jblas/operator/Tester.kt b/src/main/kotlin/com/dedztbh/kuantum/jblas/operator/Tester.kt new file mode 100644 index 0000000..122290a --- /dev/null +++ b/src/main/kotlin/com/dedztbh/kuantum/jblas/operator/Tester.kt @@ -0,0 +1,210 @@ +package com.dedztbh.kuantum.jblas.operator + +import com.dedztbh.kuantum.common.Config +import com.dedztbh.kuantum.common.readInt +import com.dedztbh.kuantum.common.upperBound +import com.dedztbh.kuantum.jblas.matrix.* +import kotlinx.coroutines.CoroutineScope +import kotlinx.coroutines.async +import kotlin.math.sqrt +import kotlin.random.Random + +/** + * Created by DEDZTBH on 2020/09/22. + * Project KuantumCircuitSim + */ + +fun getJointState(initState: String, jointStateSize: Int) = + if (initState.isBlank()) + CMatrix(jointStateSize, 1).apply { + put(0, 0, 1.0) + } + else loadCsv(initState) + +class Tester(config: Config, scope: CoroutineScope) : TFinder(config, scope) { + /** 2^N by 1 column vector */ + var jointState = getJointState(config.init_state, jointStateSize) + var hasMeasGate = false + + override suspend fun runCmd(cmd: String) = when (cmd) { + "MEASURE" -> { + val i = readInt() + val results = opMatrix * jointState + val probs = mutableListOf(0.0) + val labels = mutableListOf() + for (j in 0 until jointStateSize) { + val amag2 = results.get(0, j).magnitude2 + if (amag2 > 0) { + labels.add(j) + probs.add(amag2 + probs.last()) + } + } + println("Measurement(s):") + repeat(i) { + val r = Random.nextDouble() + val index = probs.upperBound(r) + println(allssket[labels[index]]) + } + println() + 0 + } + "MEASALL" -> { + jointState = opMatrix * jointState + opMatrix = IN2 + val probs = mutableListOf(0.0) + val labels = mutableListOf() + for (j in 0 until jointStateSize) { + val amag2 = jointState.get(j, 0).magnitude2 + if (amag2 > 0) { + labels.add(j) + probs.add(amag2 + probs.last()) + } + } + val index = probs.upperBound(Random.nextDouble()) + val collapseState = labels[index] + println("MeasAll: ${allssket[collapseState]}") + jointState = CMatrix(jointStateSize, 1).apply { + put(collapseState, 0, 1.0) + } + hasMeasGate = true + 0 + } + "MEASONE" -> { + val i = readInt() + jointState = opMatrix * jointState + opMatrix = IN2 + var prob = 0.0 + for (j in 0 until jointStateSize) { + if (allss[j][i] == '0') { + prob += jointState.get(j, 0).magnitude2 + } + } + val c = if (Random.nextDouble() < prob) '0' else { + prob = 1.0 - prob + '1' + } + println("MeasOne: Qubit #$i is |$c>") + for (j in 0 until jointStateSize) { + if (allss[j][i] != c) { + jointState.put(j, 0, 0.0, 0.0) + } + } + jointState.muli(1.0 / sqrt(prob)) + hasMeasGate = true + 0 + } + else -> super.runCmd(cmd) + } + + override suspend fun done() { + if (!hasMeasGate) super.done() + } + + override suspend fun printResult() { + if (hasMeasGate) { + println("Circuit matrix unavailable due to MeasAll/MeasOne command(s).") + } else { + super.printResult() + } + println("\nFinal state: ") + println("Init ${allssket[0]}") + (opMatrix * jointState).printFancy(allssket) + } +} + +class PTester(config: Config, scope: CoroutineScope) : PTFinder(config, scope) { + /** 2^N by 1 column vector */ + var jointState = getJointState(config.init_state, jointStateSize) + var hasMeasGate = false + + override suspend fun runCmd(cmd: String) = when (cmd) { + "MEASURE" -> { + val i = readInt() + opMatrix = reduceOps() + reversedNewOps = mutableListOf(scope.async { opMatrix }) + val results = opMatrix * jointState + val probs = mutableListOf(0.0) + val labels = mutableListOf() + for (j in 0 until jointStateSize) { + val amag2 = results.get(0, j).magnitude2 + if (amag2 > 0) { + labels.add(j) + probs.add(amag2 + probs.last()) + } + } + println("Measurement(s):") + repeat(i) { + val r = Random.nextDouble() + val index = probs.upperBound(r) + println(allssket[labels[index]]) + } + println() + 0 + } + "MEASALL" -> { + jointState = reduceOps() * jointState + reversedNewOps = mutableListOf(scope.async { IN2 }) + val probs = mutableListOf(0.0) + val labels = mutableListOf() + for (j in 0 until jointStateSize) { + val amag2 = jointState.get(j, 0).magnitude2 + if (amag2 > 0) { + labels.add(j) + probs.add(amag2 + probs.last()) + } + } + val r = Random.nextDouble() + val index = probs.upperBound(r) + val collapseState = labels[index] + println("MeasAll: ${allssket[collapseState]}") + jointState = CMatrix(jointStateSize, 1).apply { + put(collapseState, 0, 1.0) + } + hasMeasGate = true + 0 + } + "MEASONE" -> { + val i = readInt() + jointState = reduceOps() * jointState + reversedNewOps = mutableListOf(scope.async { IN2 }) + var prob = 0.0 + for (j in 0 until jointStateSize) { + if (allss[j][i] == '0') { + prob += jointState.get(j, 0).magnitude2 + } + } + val c = if (Random.nextDouble() < prob) '0' else { + prob = 1.0 - prob + '1' + } + println("MeasOne: Qubit #$i is |$c>") + for (j in 0 until jointStateSize) { + if (allss[j][i] != c) { + jointState.put(j, 0, 0.0, 0.0) + } + } + jointState.muli(1.0 / sqrt(prob)) + hasMeasGate = true + 0 + } + else -> super.runCmd(cmd) + } + + override suspend fun done() { + if (!hasMeasGate) super.done() + else { + opMatrix = reduceOps() + } + } + + override suspend fun printResult() { + if (hasMeasGate) { + println("Circuit matrix unavailable due to MeasAll/MeasOne command(s).") + } else { + super.printResult() + } + println("\nFinal state: ") + println("Init ${allssket[0]}") + (opMatrix * jointState).printFancy(allssket) + } +} \ No newline at end of file diff --git a/src/main/kotlin/matrix/CMatrixIO.kt b/src/main/kotlin/matrix/CMatrixIO.kt deleted file mode 100644 index 865ed54..0000000 --- a/src/main/kotlin/matrix/CMatrixIO.kt +++ /dev/null @@ -1,114 +0,0 @@ -package matrix - -import com.github.doyaaaaaken.kotlincsv.client.CsvReader -import com.github.doyaaaaaken.kotlincsv.client.CsvWriter -import org.ejml.UtilEjml -import org.ejml.data.MatrixSparse -import org.ejml.data.MatrixType -import org.ejml.ops.MatrixIO -import java.io.* -import java.text.DecimalFormat - -/** - * Created by DEDZTBH on 2020/09/27. - * Project KuantumCircuitSim - */ - -object CMatrixIO { - fun getMatrixType(mat: CMatrix): String { - return if (mat.type == MatrixType.UNSPECIFIED) { - mat.javaClass.simpleName - } else { - mat.type.name - } - } - - fun printTypeSize(out: PrintStream, mat: CMatrix) { - if (mat is MatrixSparse) { - val m = mat as MatrixSparse - out.println( - "Type = " + getMatrixType(mat) + " , rows = " + mat.getNumRows() + - " , cols = " + mat.getNumCols() + " , nz_length = " + m.nonZeroLength - ) - } else { - out.println("Type = " + getMatrixType(mat) + " , rows = " + mat.getNumRows() + " , cols = " + mat.getNumCols()) - } - } - - fun padSpace(builder: java.lang.StringBuilder, length: Int): String { - builder.delete(0, builder.length) - for (i in 0 until length) { - builder.append(' ') - } - return builder.toString() - } - - fun saveBin(A: CMatrix, fileName: String) { - val fileStream = FileOutputStream(fileName) - val stream = ObjectOutputStream(fileStream) - try { - stream.writeObject(A) - stream.flush() - } finally { - // clean up - fileStream.use { - stream.close() - } - } - } - - fun saveCsv(A: CMatrix, fileName: String) = - CsvWriter().writeAll(A.data.asList().windowed(A.rowStride, A.rowStride), fileName) - - @Suppress("UNCHECKED_CAST") - fun loadBin(fileName: String): T { - val fileStream = FileInputStream(fileName) - val stream = ObjectInputStream(fileStream) - val ret: T - try { - ret = stream.readObject() as T - if (stream.available() != 0) { - throw RuntimeException("File not completely read?") - } - } catch (e: ClassNotFoundException) { - throw RuntimeException(e) - } - stream.close() - return ret - } - - fun loadCsv(fileName: String): CMatrix = - CMatrix(CsvReader().open(fileName) { - readAllAsSequence().map { - DoubleArray(it.size) { i -> it[i].toDouble() } - }.toMutableList().toTypedArray() - }) - - fun CMatrix.printFancy2( - out: PrintStream = System.out, - length: Int = MatrixIO.DEFAULT_LENGTH, - allssket: List - ) { -// printTypeSize(out, this) - val format = DecimalFormat("#") - val builder = StringBuilder(length) - val cols = numCols - val c = CNum() - var i = 0 - for (y in 0 until numRows) { - for (x in 0 until cols) { - get(y, x, c) - var real = UtilEjml.fancyString(c.real, format, length, 4) - var img = UtilEjml.fancyString(c.imaginary, format, length, 4) - real += padSpace(builder, length - real.length) - img = img + "i" + padSpace(builder, length - img.length) - out.print("${allssket[i++]}: $real + $img") - if (x < numCols - 1) { - out.print(" , ") - } - } - out.println() - } - } - -} \ No newline at end of file