bind ant floats

.devcontainer/gpu/devcontainer.json

@@ -9,8 +9,8 @@
 		"--cap-add=SYS_PTRACE",
 		"--security-opt",
 		"seccomp=unconfined",
-		"--gpus",
-		"all"
+		"--device",
+		"nvidia.com/gpu=all"
 	],
 	// Set *default* container specific settings.json values on container create.
 	"customizations": {

exe/main.cpp

@@ -10,6 +10,7 @@ int main() {
     std::mt19937 mt(0);
     constexpr std::uint64_t n_qubits = 20;
     std::uniform_int_distribution<std::uint64_t> dist(0, n_qubits - 1);
+    /*
     {
         using Fp = scaluq::F16;
         scaluq::StateVector<Fp> state(n_qubits);
@@ -22,30 +23,8 @@ int main() {
         std::cout << std::chrono::duration_cast<std::chrono::milliseconds>(ed - st).count()
                   << std::endl;
     }
-    {
-        using Fp = scaluq::F32;
-        scaluq::StateVector<Fp> state(n_qubits);
-        auto st = std::chrono::system_clock::now();
-        for (int i = 0; i < 10000; i++) {
-            auto x_gate = scaluq::gate::X<Fp>(dist(mt));
-            x_gate->update_quantum_state(state);
-        }
-        auto ed = std::chrono::system_clock::now();
-        std::cout << std::chrono::duration_cast<std::chrono::milliseconds>(ed - st).count()
-                  << std::endl;
-    }
-    {
-        using Fp = scaluq::F64;
-        scaluq::StateVector<Fp> state(n_qubits);
-        auto st = std::chrono::system_clock::now();
-        for (int i = 0; i < 10000; i++) {
-            auto x_gate = scaluq::gate::X<Fp>(dist(mt));
-            x_gate->update_quantum_state(state);
-        }
-        auto ed = std::chrono::system_clock::now();
-        std::cout << std::chrono::duration_cast<std::chrono::milliseconds>(ed - st).count()
-                  << std::endl;
-    }
+    */
+    /*
     {
         using Fp = scaluq::BF16;
         scaluq::StateVector<Fp> state(n_qubits);
@@ -58,5 +37,9 @@ int main() {
         std::cout << std::chrono::duration_cast<std::chrono::milliseconds>(ed - st).count()
                   << std::endl;
     }
+    */
     Kokkos::finalize();
+
+    std::cout << "Type of _Float64: " << typeid(_Float64).name() << std::endl;
+    std::cout << "Type of double: " << typeid(double).name() << std::endl;
 }

exe/main.py

@@ -1,11 +1,40 @@
-from scaluq import *
-# from scaluq.gate import S
+import random
+import time
 
 def main():
-    state = StateVector(2)
-    swap_gate = gate.Swap(0, 1)
-    print(state)
-    print(swap_gate)
+    n_qubits = 20
+
+    from scaluq.f16 import StateVector, gate
+    state = StateVector(n_qubits)
+    st = time.time()
+    for i in range(10000):
+        gate.X(random.randint(0, n_qubits-1)).update_quantum_state(state)
+    ed = time.time()
+    print(ed-st)
+
+    from scaluq.f32 import StateVector, gate
+    state = StateVector(n_qubits)
+    st = time.time()
+    for i in range(10000):
+        gate.X(random.randint(0, n_qubits-1)).update_quantum_state(state)
+    ed = time.time()
+    print(ed-st)
+
+    from scaluq.f64 import StateVector, gate
+    state = StateVector(n_qubits)
+    st = time.time()
+    for i in range(10000):
+        gate.X(random.randint(0, n_qubits-1)).update_quantum_state(state)
+    ed = time.time()
+    print(ed-st)
+
+    from scaluq.bf16 import StateVector, gate
+    state = StateVector(n_qubits)
+    st = time.time()
+    for i in range(10000):
+        gate.X(random.randint(0, n_qubits-1)).update_quantum_state(state)
+    ed = time.time()
+    print(ed-st)
 
 if __name__ == "__main__":
     main()

include/scaluq/gate/gate_factory.hpp

@@ -242,7 +242,7 @@ inline Gate<Fp> SparseMatrix(const std::vector<std::uint64_t>& targets,
         internal::vector_to_mask(targets), internal::vector_to_mask(controls), matrix_transformed);
 }
 template <FloatingPoint Fp>
-inline Gate<Fp> Probablistic(const std::vector<Fp>& distribution,
+inline Gate<Fp> Probablistic(const std::vector<double>& distribution,
                              const std::vector<Gate<Fp>>& gate_list) {
     return internal::GateFactory::create_gate<internal::ProbablisticGateImpl<Fp>>(distribution,
                                                                                   gate_list);

include/scaluq/gate/gate_matrix.hpp

@@ -150,11 +150,11 @@ namespace internal {
 template <FloatingPoint Fp>
 void bind_gate_gate_matrix_hpp(nb::module_& m) {
     DEF_GATE(OneTargetMatrixGate, Fp, "Specific class of one-qubit dense matrix gate.")
-        .def("matrix", [](const OneTargetMatrixGate<double>& gate) { return gate->matrix(); });
+        .def("matrix", [](const OneTargetMatrixGate<Fp>& gate) { return gate->matrix(); });
     DEF_GATE(TwoTargetMatrixGate, Fp, "Specific class of two-qubit dense matrix gate.")
-        .def("matrix", [](const TwoTargetMatrixGate<double>& gate) { return gate->matrix(); });
+        .def("matrix", [](const TwoTargetMatrixGate<Fp>& gate) { return gate->matrix(); });
     DEF_GATE(SparseMatrixGate, Fp, "Specific class of sparse matrix gate.")
-        .def("matrix", [](const SparseMatrixGate<double>& gate) { return gate->get_matrix(); })
+        .def("matrix", [](const SparseMatrixGate<Fp>& gate) { return gate->get_matrix(); })
         .def("sparse_matrix",
              [](const SparseMatrixGate<Fp>& gate) { return gate->get_sparse_matrix(); });
     DEF_GATE(DenseMatrixGate, Fp, "Specific class of dense matrix gate.")

include/scaluq/gate/param_gate_factory.hpp

@@ -45,7 +45,7 @@ inline ParamGate<Fp> ParamPauliRotation(const PauliOperator<Fp>& pauli,
 }
 template <FloatingPoint Fp>
 inline ParamGate<Fp> ParamProbablistic(
-    const std::vector<Fp>& distribution,
+    const std::vector<double>& distribution,
     const std::vector<std::variant<Gate<Fp>, ParamGate<Fp>>>& gate_list) {
     return internal::ParamGateFactory::create_gate<internal::ParamProbablisticGateImpl<Fp>>(
         distribution, gate_list);
@@ -85,9 +85,9 @@ void bind_gate_param_gate_factory(nb::module_& mgate) {
               "Generate general ParamGate class instance of ParamProbablistic.");
     mgate.def(
         "ParamProbablistic",
-        [](const std::vector<std::pair<Fp, std::variant<Gate<Fp>, ParamGate<Fp>>>>&
+        [](const std::vector<std::pair<double, std::variant<Gate<Fp>, ParamGate<Fp>>>>&
                prob_gate_list) {
-            std::vector<Fp> distribution;
+            std::vector<double> distribution;
             std::vector<std::variant<Gate<Fp>, ParamGate<Fp>>> gate_list;
             distribution.reserve(prob_gate_list.size());
             gate_list.reserve(prob_gate_list.size());

pyproject.toml

@@ -1,5 +1,5 @@
 [build-system]
-requires = ["scikit-build-core >= 0.4.3", "nanobind >= 2.0.0", "typing_extensions >= 4.0.0", "setuptools >= 64", "setuptools_scm >= 8"]
+requires = ["scikit-build-core >= 0.4.3", "nanobind @ git+https://github.com/hpkfft/nanobind.git@floatingpoint", "typing_extensions >= 4.0.0", "setuptools >= 64", "setuptools_scm >= 8"]
 build-backend = "scikit_build_core.build"
 
 [project]
@@ -63,5 +63,4 @@ doc = [
     "sphinxcontrib-napoleon == 0.7",
     "sphinx_rtd_theme == 2.0.0",
     "sphinx-math-dollar == 1.2.1",
-    "nanobind == 2.0.0"
 ]

python/CMakeLists.txt

@@ -21,7 +21,20 @@ nanobind_add_stub(
     MARKER_FILE scaluq/py.typed
     VERBOSE
 )
-foreach(FLOAT IN ITEMS f64 f32)
+set(FLOATS)
+if(SCALUQ_FLOAT16)
+    list(APPEND FLOATS f16)
+endif()
+if(SCALUQ_FLOAT32)
+    list(APPEND FLOATS f32)
+endif()
+if(SCALUQ_FLOAT64)
+    list(APPEND FLOATS f64)
+endif()
+if(SCALUQ_BFLOAT16)
+    list(APPEND FLOATS bf16)
+endif()
+foreach(FLOAT IN LISTS FLOATS)
     nanobind_add_stub(
         scaluq_stub
         INSTALL_TIME

python/binding.cpp

@@ -31,17 +31,17 @@ NAMESPACE_BEGIN(NB_NAMESPACE)
 NAMESPACE_BEGIN(detail)
 
 template <typename T>
-struct type_caster<scaluq::Complex<T>> {
-    NB_TYPE_CASTER(scaluq::Complex<T>, const_name("complex"))
+struct type_caster<Complex<T>> {
+    NB_TYPE_CASTER(Complex<T>, const_name("complex"))
 
     template <bool Recursive = true>
     bool from_python(handle src, uint8_t flags, cleanup_list* cleanup) noexcept {
         (void)flags;
         (void)cleanup;
 
         if (PyComplex_Check(src.ptr())) {
-            value = scaluq::Complex<T>((T)PyComplex_RealAsDouble(src.ptr()),
-                                       (T)PyComplex_ImagAsDouble(src.ptr()));
+            value = Complex<T>((T)PyComplex_RealAsDouble(src.ptr()),
+                               (T)PyComplex_ImagAsDouble(src.ptr()));
             return true;
         }
 
@@ -59,7 +59,7 @@ struct type_caster<scaluq::Complex<T>> {
 
         make_caster<T> caster;
         if (caster.from_python(src, flags, cleanup)) {
-            value = scaluq::Complex<T>(caster.operator cast_t<T>());
+            value = Complex<T>(caster.operator cast_t<T>());
             return true;
         }
 
@@ -75,6 +75,26 @@ struct type_caster<scaluq::Complex<T>> {
     }
 };
 
+// template <>
+// struct dtype_traits<scaluq::F16> {
+//     static constexpr dlpack::dtype value{
+//         (uint8_t)dlpack::dtype_code::Float,  // type code
+//         16,                                  // size in bits
+//         1                                    // lanes (simd), usually set to 1
+//     };
+//     static constexpr auto name = const_name("float16");
+// };
+
+// template <>
+// struct dtype_traits<scaluq::BF16> {
+//     static constexpr dlpack::dtype value{
+//         (uint8_t)dlpack::dtype_code::Float,  // type code
+//         16,                                  // size in bits
+//         1                                    // lanes (simd), usually set to 1
+//     };
+//     static constexpr auto name = const_name("bfloat16");
+// };
+
 NAMESPACE_END(detail)
 NAMESPACE_END(NB_NAMESPACE)
 
@@ -116,16 +136,16 @@ NB_MODULE(scaluq_core, m) {
     internal::bind_gate_gate_hpp_without_precision(m);
 
 #ifdef SCALUQ_FLOAT16
-    bind_on_precision<std::float16_t>(m, "f16");
+    bind_on_precision<F16>(m, "f16");
 #endif
 #ifdef SCALUQ_FLOAT32
-    bind_on_precision<std::float32_t>(m, "f32");
+    bind_on_precision<F32>(m, "f32");
 #endif
 #ifdef SCALUQ_FLOAT64
-    bind_on_precision<std::float64_t>(m, "f64");
+    bind_on_precision<F64>(m, "f64");
 #endif
 #ifdef SCALUQ_BFLOAT16
-    bind_on_precision<std::bfloat16_t>(m, "bf16");
+    bind_on_precision<BF16>(m, "bf16");
 #endif
 
     initialize();

src/operator/operator.cpp

@@ -17,7 +17,7 @@ std::string Operator<Fp>::to_string() const {
 
 FLOAT(Fp)
 void Operator<Fp>::add_operator(PauliOperator<Fp>&& mpt) {
-    _is_hermitian &= mpt.coef().imag() == Fp{0.};
+    _is_hermitian &= mpt.coef().imag() == Fp{0};
     if (![&] {
             const auto& target_list = mpt.target_qubit_list();
             if (target_list.empty()) return true;

src/operator/pauli_operator.cpp

@@ -265,7 +265,8 @@ internal::ComplexMatrix<Fp> PauliOperator<Fp>::get_matrix_ignoring_coef() const
     std::uint64_t matrix_dim = 1ULL << _ptr->_pauli_id_list.size();
     internal::ComplexMatrix<Fp> mat = internal::ComplexMatrix<Fp>::Zero(matrix_dim, matrix_dim);
     for (std::uint64_t index = 0; index < matrix_dim; index++) {
-        const StdComplex<Fp> sign = Fp{1 - 2 * (Kokkos::popcount(index & phase_mask) % 2)};
+        const StdComplex<Fp> sign =
+            static_cast<Fp>(1 - 2 * (Kokkos::popcount(index & phase_mask) % 2));
         mat(index, index ^ flip_mask) = rot[rot90_count % 4] * sign;
     }
     return mat;