Add gate-batched cache-local runner across SIMD backends

apps/BUILD

@@ -32,6 +32,16 @@ cc_binary(
     ],
 )
 
+cc_binary(
+    name = "qsim_gate_batch",
+    srcs = ["qsim_gate_batch.cc"],
+    copts = gsframe_copts,
+    deps = [
+        "//lib:run_qsim_gate_batch",
+        "//lib:run_qsim_lib",
+    ],
+)
+
 cc_binary(
     name = "qsim_von_neumann",
     srcs = ["qsim_von_neumann.cc"],

apps/Makefile

@@ -1,6 +1,8 @@
 CXX_TARGETS = $(shell find . -maxdepth 1 -name '*.cc')
 CXX_TARGETS := $(CXX_TARGETS:%.cc=%.x)
 
+OUTPUT_DIR ?= .
+
 CUDA_TARGETS = $(shell find . -maxdepth 1 -name '*cuda.cu')
 CUDA_TARGETS := $(CUDA_TARGETS:%cuda.cu=%cuda.x)
 
@@ -16,6 +18,12 @@ HIP_TARGETS := $(HIP_TARGETS:%cuda.cu=%hip.x)
 .PHONY: qsim
 qsim: $(CXX_TARGETS)
 
+.PHONY: gate-batch
+gate-batch:
+	mkdir -p "$(OUTPUT_DIR)"
+	$(CXX) -o "$(OUTPUT_DIR)/qsim_gate_batch.x" qsim_gate_batch.cc $(CXXFLAGS) $(ARCHFLAGS)
+	$(CXX) -o "$(OUTPUT_DIR)/qsim_base.x" qsim_base.cc $(CXXFLAGS) $(ARCHFLAGS)
+
 .PHONY: qsim-cuda
 qsim-cuda: $(CUDA_TARGETS)
 

apps/qsim_gate_batch.cc

@@ -0,0 +1,124 @@
+// Benchmark app for the proposal-faithful gate-batch runner.
+// -f is the per-batch max_fused_size (0 = apply
+// raw gates without fusing, 2-3 = proposal's suggestion).
+
+#include <unistd.h>
+
+#ifdef _OPENMP
+# include <omp.h>
+#endif
+
+#include <limits>
+#include <string>
+
+#include "../lib/circuit_qsim_parser.h"
+#include "../lib/formux.h"
+#include "../lib/fuser_mqubit.h"
+#include "../lib/gates_qsim.h"
+#include "../lib/io_file.h"
+#include "../lib/operation.h"
+#include "../lib/run_qsim_gate_batch.h"
+#include "../lib/seqfor.h"
+#include "../lib/simmux.h"
+#include "../lib/util_cpu.h"
+
+struct Options {
+  std::string circuit_file;
+  unsigned maxtime = std::numeric_limits<unsigned>::max();
+  unsigned seed = 1;
+  unsigned num_threads = 1;
+  unsigned max_fused_size = 3;
+  unsigned block_qubits = 19;
+  unsigned verbosity = 0;
+};
+
+Options GetOptions(int argc, char* argv[]) {
+  constexpr char usage[] = "usage:\n  ./qsim_gate_batch -c circuit "
+                           "-d maxtime -s seed -t threads "
+                           "-f max_fused_size -l block_qubits -v verbosity\n";
+  Options opt;
+  int k;
+  while ((k = getopt(argc, argv, "c:d:s:t:f:l:v:")) != -1) {
+    switch (k) {
+      case 'c': opt.circuit_file = optarg; break;
+      case 'd': opt.maxtime = std::atoi(optarg); break;
+      case 's': opt.seed = std::atoi(optarg); break;
+      case 't': opt.num_threads = std::atoi(optarg); break;
+      case 'f': opt.max_fused_size = std::atoi(optarg); break;
+      case 'l': opt.block_qubits = std::atoi(optarg); break;
+      case 'v': opt.verbosity = std::atoi(optarg); break;
+      default: qsim::IO::errorf(usage); exit(1);
+    }
+  }
+  if (opt.circuit_file.empty()) {
+    qsim::IO::errorf(usage);
+    exit(1);
+  }
+  return opt;
+}
+
+template <typename StateSpace, typename QubitMappedState>
+void PrintAmplitudes(unsigned num_qubits, const StateSpace& state_space,
+                     const QubitMappedState& state) {
+  static constexpr char const* bits[8] = {
+    "000", "001", "010", "011", "100", "101", "110", "111",
+  };
+  uint64_t size = std::min(uint64_t{8}, uint64_t{1} << num_qubits);
+  unsigned s = 3 - std::min(unsigned{3}, num_qubits);
+  for (uint64_t i = 0; i < size; ++i) {
+    auto a = state.GetAmpl(state_space, i);
+    qsim::IO::messagef("%s:%16.8g%16.8g%16.8g\n",
+                       bits[i] + s, std::real(a), std::imag(a), std::norm(a));
+  }
+}
+
+int main(int argc, char* argv[]) {
+  using namespace qsim;
+
+  auto opt = GetOptions(argc, argv);
+
+#ifdef _OPENMP
+  omp_set_num_threads(opt.num_threads);
+#endif
+
+  Circuit<Operation<float>> circuit;
+  if (!CircuitQsimParser<IOFile>::FromFile(opt.maxtime, opt.circuit_file,
+                                           circuit)) {
+    return 1;
+  }
+
+  struct Factory {
+    Factory(unsigned num_threads) : num_threads(num_threads) {}
+    using Simulator = qsim::Simulator<For>;
+    using StateSpace = Simulator::StateSpace;
+    StateSpace CreateStateSpace() const { return StateSpace(num_threads); }
+    Simulator CreateSimulator() const { return Simulator(num_threads); }
+    unsigned num_threads;
+  };
+
+  using Simulator = Factory::Simulator;
+  using StateSpace = Simulator::StateSpace;
+  using State = StateSpace::State;
+  using Fuser = MultiQubitGateFuser<IO>;
+  using SeqSimulator = qsim::Simulator<SequentialFor>;
+  using Runner = QSimGateBatchRunner<IO, Fuser, Factory, SeqSimulator>;
+
+  StateSpace state_space = Factory(opt.num_threads).CreateStateSpace();
+  QubitMappedState<State> state(state_space.Create(circuit.num_qubits));
+  if (state_space.IsNull(state.state)) {
+    IO::errorf("not enough memory: is the number of qubits too large?\n");
+    return 1;
+  }
+  state_space.SetStateZero(state.state);
+
+  Runner::Parameter param;
+  param.max_fused_size = opt.max_fused_size;
+  param.block_qubits = opt.block_qubits;
+  param.verbosity = opt.verbosity;
+
+  if (Runner::Run(param, Factory(opt.num_threads), circuit, state)) {
+    PrintAmplitudes(circuit.num_qubits, state_space, state);
+  }
+
+  return 0;
+}

lib/BUILD

@@ -564,6 +564,21 @@ cc_library(
     ],
 )
 
+cc_library(
+    name = "run_qsim_gate_batch",
+    hdrs = [
+        "qubit_mapped_state.h",
+        "qubit_remap.h",
+        "run_qsim_gate_batch.h",
+    ],
+    deps = [
+        ":gate",
+        ":matrix",
+        ":operation_base",
+        ":util",
+    ],
+)
+
 cc_library(
     name = "run_qsimh",
     hdrs = ["run_qsimh.h"],

lib/qubit_mapped_state.h

@@ -0,0 +1,122 @@
+// Copyright 2026 Google LLC. All Rights Reserved.
+//
+// Licensed under the Apache License, Version 2.0 (the "License");
+// you may not use this file except in compliance with the License.
+// You may obtain a copy of the License at
+//
+//     https://www.apache.org/licenses/LICENSE-2.0
+//
+// Unless required by applicable law or agreed to in writing, software
+// distributed under the License is distributed on an "AS IS" BASIS,
+// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+// See the License for the specific language governing permissions and
+// limitations under the License.
+
+#ifndef QUBIT_MAPPED_STATE_H_
+#define QUBIT_MAPPED_STATE_H_
+
+#include <cstddef>
+#include <cstdint>
+#include <utility>
+#include <vector>
+
+#include "qubit_remap.h"
+
+namespace qsim {
+
+// Bidirectional logical<->physical qubit map; the two arrays are inverse
+// permutations of each other at all times. "Physical position p" means
+// amplitude-index bit p of the state storage.
+class QubitLayout {
+ public:
+  explicit QubitLayout(unsigned num_qubits)
+      : logical_to_physical_(num_qubits),
+        physical_to_logical_(num_qubits) {
+    for (unsigned i = 0; i < num_qubits; ++i) {
+      logical_to_physical_[i] = i;
+      physical_to_logical_[i] = i;
+    }
+  }
+
+  std::size_t NumQubits() const { return logical_to_physical_.size(); }
+
+  unsigned PhysicalPositionOf(unsigned logical_qubit) const {
+    return logical_to_physical_[logical_qubit];
+  }
+
+  unsigned LogicalQubitAt(unsigned physical_position) const {
+    return physical_to_logical_[physical_position];
+  }
+
+  const std::vector<unsigned>& LogicalToPhysical() const {
+    return logical_to_physical_;
+  }
+
+  uint64_t PhysicalIndexOf(uint64_t logical_index) const {
+    uint64_t physical_index = 0;
+    for (unsigned q = 0; q < NumQubits(); ++q) {
+      if ((logical_index >> q) & 1) {
+        physical_index |= uint64_t{1} << PhysicalPositionOf(q);
+      }
+    }
+    return physical_index;
+  }
+
+  bool IsIdentityAt(unsigned position) const {
+    return physical_to_logical_[position] == position;
+  }
+
+  void SwapPositions(unsigned position_a, unsigned position_b) {
+    const auto qubit_a = physical_to_logical_[position_a];
+    const auto qubit_b = physical_to_logical_[position_b];
+    std::swap(physical_to_logical_[position_a],
+              physical_to_logical_[position_b]);
+    std::swap(logical_to_physical_[qubit_a],
+              logical_to_physical_[qubit_b]);
+  }
+
+  // Builds one disjoint restoration pass and updates this layout to match.
+  bool BuildIdentityRestorationSwaps(std::vector<QubitSwap>& swaps) {
+    swaps.clear();
+    std::vector<char> position_touched(NumQubits(), 0);
+
+    for (unsigned position = 0; position < NumQubits(); ++position) {
+      if (IsIdentityAt(position)) continue;
+      const auto paired_position = PhysicalPositionOf(position);
+      if (position_touched[position] ||
+          position_touched[paired_position]) {
+        continue;
+      }
+
+      swaps.emplace_back(position, paired_position);
+      position_touched[position] = 1;
+      position_touched[paired_position] = 1;
+      SwapPositions(position, paired_position);
+    }
+    return !swaps.empty();
+  }
+
+ private:
+  std::vector<unsigned> logical_to_physical_;
+  std::vector<unsigned> physical_to_logical_;
+};
+
+// Owns physical state storage together with its logical-qubit interpretation.
+// A newly created mapped state starts in canonical one-to-one qubit order.
+template <typename State>
+struct QubitMappedState {
+  explicit QubitMappedState(State&& state_storage)
+      : state(std::move(state_storage)), layout(state.num_qubits()) {}
+
+  template <typename StateSpace>
+  auto GetAmpl(const StateSpace& state_space, uint64_t logical_index) const {
+    return state_space.GetAmpl(state, layout.PhysicalIndexOf(logical_index));
+  }
+
+  State state;
+  QubitLayout layout;
+};
+
+}  // namespace qsim
+
+#endif  // QUBIT_MAPPED_STATE_H_