# Scaling up HBM Efficiency of Top-K SpMV for Approximate Embedding Similarity
on FPGAs

Alberto Parravicini1, Luca Giuseppe Cellamare2, Marco Siracusa2 Marco D.
Santambrogio1 1{alberto.parravicini<EMAIL_ADDRESS>2{lucagiuseppe.cellamare<EMAIL_ADDRESS>Politecnico di
Milano, DEIB, Milan, Italy

###### Abstract

Top-K SpMV is a key component of similarity-search on sparse embeddings. This
sparse workload does not perform well on general-purpose NUMA systems that
employ traditional caching strategies. Instead, modern FPGA accelerator cards
have a few tricks up their sleeve. We introduce a Top-K SpMV FPGA design that
leverages reduced precision and a novel packet-wise CSR matrix compression,
enabling custom data layouts and delivering bandwidth efficiency often
unreachable even in architectures with higher peak bandwidth. With HBM-based
boards, we are 100x faster than a multi-threaded CPU implementation and 2x
faster than a GPU with 20% higher bandwidth, with 14.2x higher power-
efficiency.

###### Index Terms:

FPGA, SpMV, Approximate Computing, Hardware Acceleration, HBM

##  I Introduction

Information Retrieval (IR) and recommender system process an always-increasing
amount of data, often with strong real-time constraint, to suggest products,
movies, news articles to billions of users. Top-K Sparse matrix-vector
multiplication (Top-K SpMV) is a key building block of high-performance
similarity-search applications found in recommender systems that store items
or documents as sparse embeddings, short numerical representations usually
obtained through a neural network [1, 2]. Sparse embeddings guarantee low
memory footprint and reduced execution time, often without decreasing
accuracy, as they capture important information while filtering out noise [3].
Top-K SpMV computes the highest K values of the product between a sparse
matrix (a matrix where only non-zero entries are stored) and a dense vector.
In our case, it matches an input embedding against a collection of sparse
embeddings and finds the K most similar ones (Figure 1).

The computation of Top-K SpMV is extremely memory intensive and presents
sequential, indirect, and fully random memory accesses, making it unsuitable
for general-purpose architectures with traditional caching policies. Improving
its performance through a custom hardware design demands a theoretically sound
approach to optimize the data-access strategy. FPGAs can fully leverage
application-specific data representations thanks to reduced-precision
arithmetic and a customizable memory subsystem. In agreement with the Roofline
methodology [4], this approach is key to improving operational intensity and
performance in memory-starved computations. Exact numerical accuracy for large
values of K is not critical, as long as the most similar embeddings are
retrieved: compared to CPUs or GPUs, FPGAs can leverage optimized arbitrary-
precision arithmetic and offer fine-grained control over the desired target
accuracy, providing better performance, lower resource utilization, and lower
power consumption. To fully utilize off-chip bandwidth, we devise an
approximation scheme based on matrix partitioning that enables a flexible
multi-core design with a single HBM channel per core and URAM for caching
high-reuse or random access data (section III-A). On an Alveo U280 accelerator
card, our design accesses 32 HBM2 pseudo-channels with a total bandwidth of
460 GB/s. Then, we maximize bandwidth efficiency through a novel streaming
sparse matrix format that leverages coalesced HBM accesses and reduced
numerical precision, and is oblivious to the matrix non-zero entries
distribution (section III-B). We propose a novel multi-core streaming FPGA
design to approximate the computation of Top-K SpMV 111Code is available at
github.com/AlbertoParravicini/approximate-spmv-topk and present the following
contributions:

Figure 1: Top-K SpMV between a sparse matrix $\mathbf{A}$ (in our case, a
collection of sparse embeddings) and a dense vector $\mathbf{x}$ (a dense
embedding), with notation as in Section III.

  * •

A multi-core Top-K SpMV FPGA design that leverages a novel compressed matrix
representation and algorithmic approximations (section III) and to effectively
exploit the available HBM bandwidth (Section IV).

  * •

A performance evaluation of our FPGA design, showing 100x faster execution
versus a state-of-the-art CPU implementation. We are also 2x faster than a GPU
with 20% higher bandwidth, with 14.2x higher power-efficiency and no
significant loss of Top-K precision (Section V).

##  II Related Work

To the best of our knowledge, no prior work optimizes the computation of Top-K
SpMV on FPGA or GPU, although existing research covers approximation
techniques on FPGA for dense matrix multiplications and deep learning [5, 6].
On CPUs, sparse_dot_topn is a C++ multi-threaded implementation that employs
the Compressed Sparse Row (CSR) format and is specialized for sparse
embeddings and documents similarity [1]. That said, the CPU performance in
this computation is inherently held back by their limited memory bandwidth and
the inability to consistently perform quick random accesses, as there are no
guarantees that requested values have not been evicted from cache. Although no
existing work specifically optimizes Top-K SpMV on GPU, there exist multiple
highly-optimized implementations of SpMV [7, 8, 9, 10]. While modern GPUs
provide higher memory bandwidth than FPGA accelerator cards (from 549 GB/s of
the Nvidia P100 to the 1555 GB/s of the Nvidia A100), existing SpMV
implementations are often unable to utilize the available bandwidth fully
[11]. Moreover, optimizations for reduced precision data-types are currently
limited to half-precision floating-point, with no support for reduced-
precision fixed-point arithmetic [12, 13] although recent work has explored
heuristics to mix single and double precision floating-point arithmetic [14].

Existing research on SpMV optimization on FPGA focuses on sparse matrix
compression [15], moving the bottleneck from memory accesses to the input data
decompression. Umuroglu et al. [16] maximize the time in which values are kept
in a fast local cache hierarchy. The efficient use of HBM on FPGA also has
received the academic community’s attention. Sadi et al. [17] propose an SpMV
FPGA implementation that achieves significant speedup leveraging HBM and a
data-compression scheme to reduce off-chip traffic. However, directly
employing off-the-shelf SpMV FPGA designs (e.g. [18]) is undesirable: it would
waste memory to store for the entire output vector (which can contain millions
of values, Section III) and demand unnecessary computation to sort the output
afterward.

##  III Theoretical Contributions

Figure 2: Simplified approximation scheme for Top-K SpMV. No errors occur if
all partitions contain less than $\mathbf{k}$ Top-K values. TABLE I: Estimated
precision of Top-K indices for increasing number of partitions. We average
results of 1000 tests.

Number of matrix rows |  Number of partitions |  |  Top-K Value, with $\mathbf{k}=8$
---|---|---|---
|  8 |  16 |  32 |  50 |  75 |  100
$N=10^{6}$ | $\mathbf{c}=16$ |  | $1$ | $1$ | $0.999$ | $0.998$ | $0.983$ | $0.942$
$\mathbf{c}=28$ |  | $1$ | $1$ | $1$ | $0.999$ | $0.999$ | $0.996$
$\mathbf{c}=32$ |  | $1$ | $1$ | $1$ | $0.999$ | $0.999$ | $0.997$
$N=10^{7}$ | $\mathbf{c}=16$ |  | $1$ | $1$ | $1$ | $0.999$ | $0.986$ | $0.947$
$\mathbf{c}=28$ |  | $1$ | $1$ | $1$ | $0.999$ | $0.999$ | $0.995$
$\mathbf{c}=32$ |  | $1$ | $1$ | $1$ | $0.999$ | $0.998$ | $0.998$

Two key theoretical contributions enable our FPGA hardware design. First, we
introduce a partitioning scheme that splits and approximates Top-K SpMV on
independent cores, giving access to the entire HBM bandwidth (Section III-A).
Then, we present Block-Streaming CSR (BS-CSR), a sparse matrix format
optimized for streaming coalesced memory transactions and reduced-precision
arithmetic (Section III-B).

Given a sparse matrix $\mathbf{A}$ of size $N\times M$, and a dense vector
$\mathbf{x}$ of size $M\times 1$, the result $\mathbf{y}$ of the SpMV
$\mathbf{y}=\mathbf{A}\mathbf{x}$ is a vector $N\times 1$. Top-K SpMV does not
track the entire $\mathbf{y}$, but only its highest K values (and the
corresponding indices). If $\mathbf{A}$ and $\mathbf{x}$ are $L2$-normalized
embeddings, Top-K SpMV retrieves the rows of $\mathbf{A}$ with highest cosine
similarity to $\mathbf{x}$. In our task, $N$ is the size of the embedding
collection (usually millions of entries), while $M$, the dense embedding size,
has only a few hundred values [19]. $\mathbf{x}$ is stored in URAM to perform
random accesses in a single clock cycle, which is infeasible with an
unconstrained $M$. Similarly, we track just the Top-K entries of $\mathbf{y}$
on a LUT scratchpad instead of performing data-dependent write-backs on HBM
and sharing bandwidth with read-channels, which would undermine the maximum
memory throughput. $\mathbf{x}$ is dense as many sparsification algorithms
operate on dense matrices (e.g. batches of $\mathbf{A}$) and cannot
efficiently sparsify a single vector [20, 21].

###  III-A Top-K SPMV Approximation

In Top-K similarity-search applications, it is often acceptable to trade
maximum accuracy on large values of $K$ for lower execution time and power
consumption. Instead of computing the total Top-K values on the matrix A, we
partition it in $\mathbf{c}$ sub-matrices with $N/\mathbf{c}$ rows each. Each
sub-matrix is processed by an independent FPGA core, which computes the top
$\mathbf{k}<K$ results for the partition (with $\mathbf{k}\cdot\mathbf{c}\geq
K$). We obtain $\mathbf{k}\cdot\mathbf{c}$ results, which represents an
approximation of the original Top-K (Figure 2). As $\mathbf{c}$ increases, so
does the approximation accuracy. As we always retrieve the top $\mathbf{k}$
values, the approximation does not affect the best-ranked rows. Equation 1
presents a closed-form expression of the expected value of correctly retrieved
Top-K values (i.e. the precision); intuitively, Equation 1 compares the number
of permutations with more than $\mathbf{k}$ Top-K values in a partition with
the total number of Top-K permutations. We estimate $\mathbb{E}[P]$ through a
Monte Carlo simulation for different $K$ and number of partitions, with $K$
from 8 to 100, common thresholds in IR [22]. Having at least 16 partitions
guarantees a minimal loss of precision, even for large matrices.

$\mathbb{E}[P]\approx\frac{1}{K-1}\sum_{K_{i}=1}^{K}{1-\frac{\mathbf{c}\cdot\sum_{\mathbf{k}_{i}=\mathbf{k}+1}^{min(K_{i},\lfloor
N/\mathbf{c}\rfloor)}{\binom{\lfloor
N/\mathbf{c}\rfloor}{\mathbf{k}_{i}}}}{\binom{N}{K_{i}}}}$ (1)

###  III-B The Block-Streaming CSR Matrix Layout

The matrix $\mathbf{A}$ is stored in a sparse format to save only non-zero
values and lower its memory footprint. Different storage techniques are
possible, depending on the non-zero elements distribution and the desired type
of access patterns.

The common CSR format [23] is unsuitable for fully-pipelined streaming FPGA
designs as it contains data dependencies to access matrix values, and hiding
the memory controller latency is not trivial, given the lack of a built-in
hardware prefetcher. Instead, the Coordinate (COO) layout uses three equally
sized arrays to store, for each matrix entry, its two coordinates and the
entry itself’s value. COO allows streaming iterations with burst memory
transactions on non-zero matrix entries to saturate bandwidth, as the
architecture does not have to perform data-dependent memory accesses
determined by the number of non-zero values of each row, as in the CSR format.
COO allows coalesced memory transactions of multiple non-zero entries, but
requires redundant storage of coordinate values, limiting the overall
operational intensity.

Figure 3: Naïve COO allows only 5 non-zero entries in a 512-bits packet. Using
fewer bits for y and val (assuming y $<1024$ and val stored as 20-bits fixed-
point) permits only 3 additional non-zeros per packet. BS-CSR drastically
reduces the space required for ptr, enabling 15 non-zeros per packet.

We propose BS-CSR, a new sparse matrix storage format that combines the low
memory footprint of CSR with the streaming properties of COO, and is enabled
by reduced precision data-types. As the architecture of the Alveo U280 HBM
memory controllers incentivizes memory transactions of data packets (or
blocks) between 256 and 512 bits [24], we build each packet in BS-CSR as an
independent CSR partition (Figure 3): values of idx and val are identical to
standard CSR, using less than 32 bits for each value if possible (in our case,
10 bits are enough for idx entries). The ptr vector tracks the cumulative non-
zero entries with respect to the packet itself: each entry of ptr is just
$\lfloor log_{2}(B)\rfloor$ bits for a packet with $B$ non-zero entries (4
bits for a packet $B\\!=\\!15$), instead of 32 bits like in a naïve CSR or COO
(as the number of rows is not bound). The $new\\_row$ bit tracks if the first
row of a packet continues the last row of the previous packet. Missing rows
are handled with placeholder $0$ values, although rows are never fully empty
in our application domain. BS-CSR does not provide an explicit reference to
row indices: this is not a limitation in streaming algorithms such as SpMV and
its variations, as we process the entire matrix sequentially and track the
current rows by counting non-zero entries in each ptr packet, plus the
$new\\_row$ bit. We can fit 2 to 3 times as many non-zero entries in each
512-bits packet and improve the operational intensity by the same amount.
Although significant work has been made to optimize SpMV on matrices with
irregular row density distribution [17], we avoid the problem entirely through
a fully-streaming data format that is oblivious to the number of non-zero
matrix entries per row. In spite of its simplicity, BS-CSR can be a life-saver
in streaming memory-bound computations (Section V-C).

##  IV The Proposed FPGA Hardware Design

Our Top-K SpMV FPGA design offers multiple low-profile cores that operate
independently. Each core processes a partition of the input matrix, as
explained in Section III-A. To guarantee maximum flexibility in terms of the
number of cores placed on the FPGA, each core uses a single HBM channel to
read the input matrix and store its output at the very end of the computation.
As cores read 512 bits per clock cycle (at 225 MHz) from their respective HBM
memory port and memory transactions happen in continuous maximum length AXI4
bursts (256 beats), our FPGA design can theoretically operate at the maximum
bandwidth offered by HBM, by coalescing transactions with BS-CSR packets
containing multiple non-zero entries, and without expensive distributed memory
controllers. Each core reads $\mathbf{A}$ in packets of size $B$. $B$ ranges
from 7 to 15, depending on the desired level of numerical precision and the
embedding size: at worst we use 32 bits for idx and val, but realistic size
bounds (e.g. idx $<1024$) allow much greater coalescing and operational
intensity (Figure 6).

Figure 4: Block diagram of our multi-core Top-K SpMV, simplified for BS-CSR
with 5 non-zero y and val entries per packet.

###  IV-A Leveraging URAM for Fast Random Access

The input vector $\mathbf{x}$ is stored in URAM. Our core performs $B$ random
accesses per cycle on $\mathbf{x}$. As each URAM bank only has 2 read ports,
we replicate $\mathbf{x}$ $\lceil B/2\rceil$ times to allow all random
accesses. This replication does not constitute a limitation of our design:
$\mathbf{x}$ represents a dense embedding whose size does not go beyond a few
hundred values in real applications. With our implementation, $\mathbf{x}$ can
have size up to 80000 (assuming a worst-case scenario with 32 bits values, 32
cores, and 8 replicas of $\mathbf{x}$ per core), given a URAM size of around
90MB.

###  IV-B Top-K SpMV Algorithm Design

The main Top-K SpMV algorithm (alg. 1) is implemented as a 4-stage pipelined
data-flow computation, decoupling memory transfer and computation. Each stage
processes the input matrix as packets of $B$ elements. The first stage
retrieves $B$ entries from $\mathbf{A}$ per cycle from HBM and URAM and
computes point-wise products. Values corresponding to the same matrix row are
aggregated in the second stage. The third stage performs book-keeping for rows
that span more than one packet and finds which rows have been fully processed
in the current packet. The final stage updates the current Top-K values for
each complete row if the row output value is in the Top-K. Our data-flow
design completely hides the Top-K update cost, which would otherwise be an
expensive part of the computation.

Each core retrieves the top $\mathbf{k}=8$ values of each matrix partition.
Higher $\mathbf{k}$ results in lower clock speed due to RAW data-dependencies
in the argmin computation, while lower $\mathbf{k}$ decreases accuracy. It is
unlikely that many rows are entirely contained in a single packet $B$: we
avoid consumption of unnecessary logic by tracking results for at most
$\mathbf{r}$ rows per packet. In our experiments, using $B/4<\mathbf{r}<B/2$
provided resource savings up to 50% with no accuracy loss.

###  IV-C Lower Precision, More Cores, Better Performance

Our hardware design processes $\mathbf{c}\cdot\\!B$ non-zeros per clock cycle.
Maximizing this performance equation is not trivial, as increasing $B$
increments FPGA resource consumption and possibly prevents the placement of 32
cores, a requirement to achieve full HBM bandwidth utilization. Indeed,
$\mathbf{c}$ is a non-linear function of $\mathbf{k}$, $\mathbf{r}$ and $B$,
with $B$ depending from $M$ and from the number of bits $V$ used for entries
of val, as in $B\cdot(\lceil log_{2}{B}\rceil+\lceil
log_{2}{M}\rceil+V)+1\\!=\\!512$. $M$ can be safely assumed in the order of
thousands. Choosing $V$, $\mathbf{k}$ and $\mathbf{r}$ requires a performance-
accuracy trade-off: based on theoretical bounds in Table I and the experiments
in Section V-D, we observe that $\mathbf{k}=8$ and $V=20$ are sufficient to
provide great accuracy. Combined with $\mathbf{r}$ between 4 and 8 we
guarantee a 32-cores design that can exploit all HBM channels without
significant accuracy loss. Table II summarizes characteristics and resource
usage of the designs in our evaluation.

Algorithm 1 Approximate BS-CSR Top-K SpMV

1:function Top-K-SpMV($bscsr\\_matrix,\mathbf{vec}$)

2: $x_{s\\_old}\leftarrow 0;\ last\\_packet\\_output\leftarrow 0$

3: for $i\leftarrow 0..\text{NNZ}/B$ do $\triangleright$ For each packet of
the matrix

4: $\triangleright$ 1\. Process BS-CSR in packets of size $B$

5: for $j\leftarrow 0..B$ do $\triangleright$ All loops of size $B$ are
Unrolled

6: $bscsr\\_packet\leftarrow bscsr\\_matrix[i]$

7: $x_{loc}[j]\leftarrow bscsr\\_packet.x[j]$

8: $val_{loc}[j]\leftarrow bscsr\\_packet.val[j]$

9: $res_{tmp}[j]\leftarrow
val_{loc}[j]\cdot\mathbf{vec}[j][bscsr\\_packet.y[j]]$

10: $\triangleright$ 2\. Aggregate partial $res_{loc}$ values

11: for $b1\leftarrow 0..B$ do $\triangleright$ $x_{loc}[-1]=0$ by convention

12: $row_{curr}\mathrel{+}=(x_{loc}[b1]\neq x_{loc}[b1-1])$

13: for $b2\leftarrow x_{loc}[b1-1]..x_{loc}[b1]$ do

14: $res_{agg}[b1]\mathrel{+}=res_{tmp}[b2]$

15: $row_{curr}\mathrel{+}=bscsr\\_packet.new\\_row-1$

16: $\triangleright$ 3\. Check if first value was split among packets

17: for $j\leftarrow 1..B$ do $\triangleright$ Find finished rows

18: $rows_{fin}[j]\leftarrow x_{loc}[j-1]\neq 0$

19: if $bscsr\\_packet.new\\_row$ then

20: $res_{agg}[0]\leftarrow last\\_packet\\_output$

21: $rows_{fin}[0]\leftarrow true$

22: else$\triangleright$ This packet continues a previous row

23: $res_{agg}[1]\mathrel{+}=last\\_packet\\_output$

24: $res_{agg}[0]\leftarrow 0;\ rows_{fin}[0]\leftarrow false$

25: $\triangleright$ 4\. Update Top-K values

26: for $j\leftarrow 0..B$ do $\triangleright$ Process only finished rows

27: if $res_{agg}[j]\geq worst_{curr}[j]\ \&\&\ rows_{fin}[j]$ then

28: $res_{loc}[j]\leftarrow res_{agg}[j]$

29: $res_{idx}[j]\leftarrow row_{curr}+j-1$

30: $argmin(res_{loc},\ worst_{curr}[j],\ worst_{idx}[j])$

31: $reset(res_{agg})$

##  V Experimental Evaluation

We employ a Xilinx Alveo U280 Accelerator Card equipped with 8 GB of HBM2
memory (460 GB/s of total bandwidth over 32 channels) and a
xcu280-fsvh2892-2L-e FPGA whose available resources are reported in Table II.
Our CPU baseline is sparse_dot_topn, a multi-threaded C++ implementation of
Top-K SpMV [1], running on two Intel Xeon Gold 6248 and 384 GB of DRAM. We
also compare our design against GPUs, using a Tesla P100 (549 GB/s of HBM
bandwidth): although we are not aware of any GPU implementation of Top-K SpMV,
we can combine a fast GPU implementation of SpMV (such as cuSPARSE [7]) with
the GPU radix sort of Thrust [25], using both single and half-precision
floating-point arithmetic. To provide a worst-case comparison, we also assume
a zero-cost GPU sorting, as if cuSPARSE already retrieved Top-K values at no
cost.

We analyze 4 FPGA design: 32 bits unsigned fixed-point (Q1.31), 25 bits
(Q1.24), 20 bits (Q1.19), and a 32-bit floating-point version (F32) (Table
II). The number of HBM channels limits the maximum number of cores to 32,
although we could easily place more cores given our design’s low resource
footprint. The CPU baseline uses floating-point arithmetic: our CPU does not
support arbitrary reduced-precision, and simulated reduced-precision fixed-
point resulted in lower performance. The experimental setup comprises 19
synthetic and real sparse embedding matrices (Table III), with different sizes
and distributions. Sizes are reported using BS-CSR as in Figure 3: if stored
as a naïve COO, they would take 3 times as much space. We simulate different
non-zero distributions (uniform and left-skewed $\Gamma$), with 20 or 40
average non-zero entries per row (a sparsity factor of 2-8%). The use of
synthetic matrices provides full control over the desired data-distribution,
as the non-zero distribution in real embedding matrices is influenced by the
chosen sparsification procedure. To the best of our knowledge, no public data-
set of sparse embedding with sizes comparable to ours (millions of rows) is
available. Instead, we sparsify the GloVe [26] embedding corpus with the
technique in [21]. Our performance is mostly unaffected by the underlying
data-distribution (Section V-A): as our FPGA design processes the matrix in a
streaming fashion, density does not affect performance. We perform each test
30 times, with different random vertices $\mathbf{x}$.

TABLE II: Resource usage, clock frequency and power consumption of our
architecture. Best values in bold.

Bit-width |  Cores |  LUT |  FF |  BRAM |  URAM |  DSP |  Clock (MHz) |  Power Cons.
---|---|---|---|---|---|---|---|---
20 bits | 32 | 38% | 35% | 20% | 33% | 7% | 253 | 34 W
25 bits | 32 | 38% | 36% | 20% | 30% | 11% | 240 | 35 W
32 bits | 32 | 35% | 33% | 20% | 27% | 17% | 249 | 35 W
32 bits, float | 32 | 44% | 37% | 20% | 26% | 19% | 204 | 45 W
Available |  | 1097419 | 2180971 | 1812 | 960 | 9020 |  |

TABLE III: Matrices in the evaluation, with $M=512$ and $1024$, and with
memory occupation using BS-CSR as in Figure 3. Distribution controls the
number of non-zeros per row.

Distribution |  Rows |  Non-zeros (min-max) |  Size (min-max, GB)
---|---|---|---
Uniform | $0.5\cdot 10^{7}$ | $10^{8}-2\cdot 10^{8}$ | $$0.4\text{\,}\mathrm{GB}$-$0.8\text{\,}\mathrm{GB}$$
$1.0\cdot 10^{7}$ | $2\cdot 10^{8}-4\cdot 10^{8}$ | $$0.8\text{\,}\mathrm{GB}$-$1.7\text{\,}\mathrm{GB}$$
$1.5\cdot 10^{7}$ | $3\cdot 10^{8}-6\cdot 10^{8}$ | $$1.2\text{\,}\mathrm{GB}$-$2.5\text{\,}\mathrm{GB}$$
$\Gamma(k=3,\theta=4/3)$ | $0.5\cdot 10^{7}$ | $9.7\cdot 10^{7}-1.97\cdot 10^{8}$ | $$0.4\text{\,}\mathrm{GB}$-$0.8\text{\,}\mathrm{GB}$$
$1.0\cdot 10^{7}$ | $1.9\cdot 10^{8}-3.95\cdot 10^{8}$ | $$0.8\text{\,}\mathrm{GB}$-$1.7\text{\,}\mathrm{GB}$$
$1.5\cdot 10^{7}$ | $2.9\cdot 10^{8}-5.92\cdot 10^{8}$ | $$1.2\text{\,}\mathrm{GB}$-$2.5\text{\,}\mathrm{GB}$$
Sparsified GloVe | $0.2\cdot 10^{7}$ | $2.4\cdot 10^{7}-4.6\cdot 10^{7}$ | $$0.1\text{\,}\mathrm{GB}$-$0.3\text{\,}\mathrm{GB}$$

Figure 5: Execution time speedup on Top-K SpMV of our FPGA design and GPU vs.
CPU. We can process over 57 billion non-zeros per second, 2x more than the GPU
SpMV.
Figure 6: Roofline model of our Top-K SpMV architecture. a Operational
intensity increase with BS-CSR. b Comparison of FPGA vs. CPU and GPU. Thanks
to BS-CSR, we provide the highest operational intensity and the best
performance.

###  V-A Execution Time

For each design and matrix, we measure the speedup against CPU and GPU
implementations and report results in Figure 5, for $K\\!=\\!100$. While both
GPU and FPGA implementations are significantly faster than the CPU baseline,
our FPGA design achieves a 2x speedup over an idealized GPU Top-K SpMV
implementation with zero-cost sorting, despite the 20% lower peak memory
bandwidth. When accounting for sorting costs on the GPU, this speedup can be
as large as 7x. We expect to provide competitive performance even against a
GPU with significantly higher memory bandwidth (e.g. Tesla A100, 1.5TB/s), and
always provide greater power-efficiency (Section V-B). Fixed-point arithmetic
guarantees higher speedups than floating-point, thanks to the pipelines’ lower
initiation interval. Reduced-precision enables packing more non-zeros per
transfer (higher $B$), resulting in increased operational intensity and better
performance. Our 32-cores design can find the Top-K values of a matrix with
$10^{7}$ rows and 200 million non-zero entries in less than 4 ms and is
suitable for real-time applications.

###  V-B Power Efficiency

Our FPGA hardware design consumes about 35W during execution, plus 40W for the
host server, measured with an external power meter monitor. The CPU
implementation consumes around 300W during execution, while the GPU requires
250W (plus 40W for the host). Our fixed-point FPGA design provides 400x higher
Performance/Watt ratio than the CPU, and 14.2x compared to the idealized GPU
implementation (7.7x when accounting for an equal host machine): we provide
higher performance without any sacrifice of power efficiency.

###  V-C Roofline Model Analysis

From the Roofline Model in Figure 6, built following the methodology in [4],
we observe how the performance (expressed as non-zeros per second) of our FPGA
design scales linearly to the total bandwidth of the available HBM channels:
this result is significant as it entails predictable performance when
deploying our design on an FPGA board with fewer (or more) HBM channels. Most
importantly, BS-CSR increases the operational intensity up to 3x compared to a
naïve COO ($B\\!=\\!15$ vs. $B\\!=\\!5$), which immediately translates to an
equivalent performance improvement. When compared against CPU and GPU, our
FPGA design is the best in terms of operational intensity and bandwidth usage,
both in absolute and percentage terms. Despite having access to 20% less
bandwidth than the GPU and considering a Top-K SpMV GPU implementation with
zero-cost sorting, we still provide the highest performance, thanks to the
increased operational intensity of BS-CSR.

Figure 7: Top-K SpMV accuracy (higher is better) for different architectures
and types of reduced-precision arithmetic.

###  V-D Approximation Accuracy Analysis

We compare the accuracy of our approximated Top-K SpMV with the exact CPU
results and with the approximated results of the GPU with half-precision
floating-point. We look at values of K from 8 to 100 and evaluate accuracy
using common Recommender System evaluation metrics such as Precision,
Kendall’s $\tau$, and NDCG [27]. Precision does not penalize out-of-order
results; the other two metrics do. Our results are in line with the
theoretical estimations in Table I, showing only a minor accuracy dip for
large K. Even 20-bit fixed-point designs provide excellent accuracy across the
board, with Precision above 97% (Figure 7). Moreover, 32-bits fixed-point
designs provide accuracy above the half-precision floating-point GPU
implementation, despite our algorithmic approximation.

##  VI Conclusion and future work

Top-K SpMV is a cornerstone of IR and recommender systems based on sparse
embeddings similarity and must be computed guaranteeing real-time latencies
and power efficiency. As these constraints are hardly met on general-purpose
architectures, we presented a novel approximate FPGA multi-core design for
Top-K SpMV that leverages HBM, fixed-point arithmetic, and our new BS-CSR
matrix layout to achieve high-performance and modularity over the FPGA
resources. Our 32-cores 20-bit FPGA design achieves 100x the performance of a
state-of-the-art CPU baseline and 2x the GPU performance with 20% higher
bandwidth. We deliver real-time results on sparse matrices with hundreds of
millions of values while keeping 14.2x higher power efficiency than the GPU.

Future work will focus on adaptive compressed matrix representations by
reconfiguring the FPGA in terms of numerical precision to guarantee desired
targets of accuracy or performance. We will also apply our design to smaller
FPGA accelerator cards: with similar memory bandwidth, the computation can be
cheaper and even more power-efficient, with no performance loss. Finally, we
will study how to apply BS-CSR to other computations and possibly other
architectures.

## References

  * [1] 2017\. [Online]. Available: https://github.com/ing-bank/sparse_dot_topn
  * [2] M. Aumüller, E. Bernhardsson, and A. Faithfull, “Ann-benchmarks: A benchmarking tool for approximate nearest neighbor algorithms,” 2018.
  * [3] G. Berend, “Sparse coding of neural word embeddings for multilingual sequence labeling,” _Transactions of the Association for Computational Linguistics_ , vol. 5, pp. 247–261, 2017.
  * [4] M. Siracusa, M. Rabozzi, E. Del Sozzo, L. Di Tucci, S. Williams, and M. D. Santambrogio, “A cad-based methodology to optimize hls code via the roofline model,” in _2020 IEEE/ACM ICCAD_ , 2020.
  * [5] H. Waris, C. Wang, W. Liu, and F. Lombardi, “Axsa: On the design of high-performance and power-efficient approximate systolic arrays for matrix multiplication,” _Journal of Signal Processing Systems_ , 8 2020.
  * [6] E. Adams, S. Venkatachalam, and S. Ko, “Energy-efficient approximate mac unit,” in _2019 IEEE ISCAS_ , 2019, pp. 1–4.
  * [7] M. Naumov, L. Chien, P. Vandermersch, and U. Kapasi, “Cusparse library,” in _GPU Technology Conference_ , 2010.
  * [8] Y. Liu and B. Schmidt, “Lightspmv: Faster csr-based sparse matrix-vector multiplication on cuda-enabled gpus,” in _2015 IEEE 26th International Conference on Application-specific Systems, Architectures and Processors (ASAP)_. IEEE, 2015, pp. 82–89.
  * [9] M. Steinberger, A. Derlery, R. Zayer, and H.-P. Seidel, “How naive is naive spmv on the gpu?” in _2016 IEEE High Performance Extreme Computing Conference (HPEC)_. IEEE, 2016, pp. 1–8.
  * [10] W. Liu and B. Vinter, “Csr5: An efficient storage format for cross-platform sparse matrix-vector multiplication,” in _Proceedings of the 29th ACM on International Conference on Supercomputing_ , 2015.
  * [11] T. Nguyen, S. Williams, M. Siracusa, C. MacLean, and N. W. D. Doerfler, “The performance and energy efficiency potential of fpgas in scientific computing,” in _IEEE/ACM PMBS_ , 2020.
  * [12] A. Abdelfattah, H. Anzt, E. G. Boman, E. Carson, T. Cojean, J. Dongarra, M. Gates, T. Grützmacher, N. J. Higham, S. Li _et al._ , “A survey of numerical methods utilizing mixed precision arithmetic,” _arXiv preprint arXiv:2007.06674_ , 2020.
  * [13] A. Haidar, S. Tomov, J. Dongarra, and N. J. Higham, “Harnessing gpu tensor cores for fast fp16 arithmetic to speed up mixed-precision iterative refinement solvers,” in _SC18_. IEEE, 2018.
  * [14] K. Ahmad, H. Sundar, and M. Hall, “Data-driven mixed precision sparse matrix vector multiplication for gpus,” _ACM Transactions on Architecture and Code Optimization (TACO)_ , vol. 16, no. 4, pp. 1–24, 2019.
  * [15] P. Grigoras, P. Burovskiy, E. Hung, and W. Luk, “Accelerating spmv on fpgas by compressing nonzero values,” in _2015 IEEE 23rd Annual International Symposium on Field-Programmable Custom Computing Machines_. IEEE, 2015, pp. 64–67.
  * [16] Y. Umuroglu and M. Jahre, “A vector caching scheme for streaming fpga spmv accelerators,” in _International Symposium on Applied Reconfigurable Computing_. Springer, 2015, pp. 15–26.
  * [17] F. Sadi, J. Sweeney, T. M. Low, J. C. Hoe, L. Pileggi, and F. Franchetti, “Efficient spmv operation for large and highly sparse matrices using scalable multi-way merge parallelization,” in _Proceedings of the 52nd Annual IEEE/ACM International Symposium on Microarchitecture_. New York, NY, USA: Association for Computing Machinery, 2019, p. 347–358.
  * [18] A. Parravicini, F. Sgherzi, and M. D. Santambrogio, “A reduced-precision streaming spmv architecture for personalized pagerank on fpga,” in _Proceedings of the 26th Asia and South Pacific Design Automation Conference_ , 2021, pp. 378–383.
  * [19] Z. Yin and Y. Shen, “On the dimensionality of word embedding,” in _Advances in Neural Information Processing Systems_ , 2018, pp. 887–898.
  * [20] J. A. Tropp and A. C. Gilbert, “Signal recovery from random measurements via orthogonal matching pursuit,” _IEEE Transactions on information theory_ , vol. 53, no. 12, pp. 4655–4666, 2007.
  * [21] J. Mairal, F. Bach, J. Ponce, and G. Sapiro, “Online dictionary learning for sparse coding,” in _Proceedings of the 26th annual international conference on machine learning_ , 2009, pp. 689–696.
  * [22] M. F. Dacrema, P. Cremonesi, and D. Jannach, “Are we really making much progress? a worrying analysis of recent neural recommendation approaches,” in _Proceedings of the 13th ACM Conference on Recommender Systems_ , 2019, pp. 101–109.
  * [23] C. Yang, A. Buluç, and J. D. Owens, “Design principles for sparse matrix multiplication on the gpu,” in _European Conference on Parallel Processing_. Springer, 2018, pp. 672–687.
  * [24] Z. Wang, H. Huang, J. Zhang, and G. Alonso, “Shuhai: Benchmarking high bandwidth memory on fpgas,” in _2020 IEEE 28th Annual International Symposium on Field-Programmable Custom Computing Machines (FCCM)_. IEEE, 2020, pp. 111–119.
  * [25] N. Bell and J. Hoberock, “Thrust: A productivity-oriented library for cuda,” in _GPU computing gems Jade edition_. Elsevier, 2012.
  * [26] J. Pennington, R. Socher, and C. D. Manning, “Glove: Global vectors for word representation,” in _Empirical Methods in Natural Language Processing (EMNLP)_ , 2014, pp. 1532–1543.
  * [27] G. Shani and A. Gunawardana, “Evaluating recommendation systems,” in _Recommender systems handbook_. Springer, 2011, pp. 257–297.


# On BGS-conjecture, CLT and measurement

Alexey A. Kryukov  Department of Mathematics & Natural Sciences, University of
Wisconsin-Milwaukee, USA<EMAIL_ADDRESS>
###### Abstract

The BGS conjecture states that the Hamiltonian of a microscopic analogue of a
classical chaotic system can be modeled by a random matrix from a Gaussian
ensemble. Here, this conjecture is considered in the context of a recently
discovered geometric relationship between classical and quantum mechanics.
Based on BGS, we conjecture that the Hamiltonian of a system whose classical
counterpart performs a random walk can be modeled by a family of independent
random matrices from the Gaussian unitary ensemble. By accepting this
conjecture, we find a relationship between the process of observation in
classical and quantum physics, and describe the boundary between the micro and
macro worlds.

##  1 Main

The linear nature of quantum mechanics poses a persistent problem for
reconciliation of the classical and quantum mechanics. The superposition
principle is foreign to classical physics. We don’t see macroscopic objects in
two different places or the cat being alive and dead at the same time.
However, such states are commonplace in the microworld. The notable attempts
to resolve the situation include non-linear modifications of the Schrödinger
equation to account for the transition to states observed in a measurement,
the De Broglie-Bohm theory of classical particles lead by a pilot-wave, an
appeal to many coexisting worlds representing the components of a
superposition and the decoherence program aiming to explain how superpositions
decohere to probabilistic mixtures. None of these attempts is generally
accepted as successful in resolving the problem. There is no experimental
evidence for needing to modify the Schrödinger equation. The simultaneous
presence of quantum and classical trajectories in the pilot-wave theory seems
redundant. The many-world approach fails to explain the world that is real and
unique to us. The decoherence program derives the laws of probability valid
for macroscopic bodies but fails to account for a specific outcome of a
measurement.

On the other hand, the Newtonian and Schrödinger dynamics have a simple
relationship that does not seem to be known nor used in the above-mentioned
attempts. Namely, the Newtonian dynamics can be identified with a constrained
Schrödinger dynamics. The latter is similar to the Newtonian dynamics of a
constrained system, e.g., a bead on a wire. However, since the Schrödinger
equation describes the dynamics in a Hilbert space of states, the constraint
is imposed on the state of the system. Consider the subset $M^{\sigma}_{3,3}$
of the Hilbert space $L_{2}(\mathbb{R}^{3})$, formed by the states

$\varphi({\bf x})=g_{{\bf a},\sigma}({\bf x})e^{i{\bf p}{\bf x}/\hbar},$ (1)

where

$g_{{\bf
a},\sigma}=\left(\frac{1}{2\pi\sigma^{2}}\right)^{3/4}e^{-\tfrac{({\bf x}-{\bf
a})^{2}}{4\sigma^{2}}}$ (2)

is the Gaussian function of a sufficiently small variance $2\sigma^{2}$
centered at a point ${\bf a}$ in the Euclidean space $\mathbb{R}^{3}$ and
${\bf p}$ is a fixed vector in $\mathbb{R}^{3}$. The state $\varphi({\bf x})$
represents a narrow wave packet with group velocity ${\bf p}/m$, where $m$ is
the mass of the particle. Let’s identify the set of all pairs $({\bf a},{\bf
p})$ with the classical phase space $\mathbb{R}^{3}\times\mathbb{R}^{3}$ of
possible positions ${\bf a}$ and momenta ${\bf p}$ of a particle. The map
$\Omega:({\bf a},{\bf p})\longrightarrow g_{{\bf a},\sigma}e^{i{\bf p}{\bf
x}/\hbar}$ identifies then the classical phase space with the submanifold
$M^{\sigma}_{3,3}$ of $L_{2}(\mathbb{R}^{3})$. The equivalence classes of
states in $L_{2}(\mathbb{R}^{3})$ differing only by a constant phase factor
$e^{i\alpha}$ form the projective space $CP^{L_{2}}$. Under the equivalence
relation, the embedded manifold $M^{\sigma}_{3,3}$ becomes a submanifold of
$CP^{L_{2}}$, denoted here by the same symbol.

The six-dimensional manifold $M^{\sigma}_{3,3}$ is embedded into the space of
states in a very special way. There are no vectors in the Hilbert space
$L_{2}(\mathbb{R}^{3})$ orthogonal to all of $M^{\sigma}_{3,3}$. Instead, the
points of $M^{\sigma}_{3,3}$ represent an overcomplete basis in the Hilbert
space [4]. Furthermore, the projective space $CP^{L_{2}}$ possesses the
Fubini-Study metric, induced by the embedding of $CP^{L_{2}}$ into the sphere
$S^{L_{2}}$ of unit-normalized states in $L_{2}(\mathbb{R}^{3})$ furnished
itself with the induced round Riemannian metric. For any two vectors
$\xi,\eta$ in $L_{2}(\mathbb{R}^{3})$ tangent to the sphere $S^{L_{2}}$ and
the corresponding vectors $X=(\mathrm{Re}\xi,\mathrm{Im}\xi)$ and
$Y=(\mathrm{Re}\eta,\mathrm{Im}\eta)$ in the realization
$L_{2R}(\mathbb{R}^{3})$ of $L_{2}(\mathbb{R}^{3})$, the Riemannian metric
$G_{\varphi}$ on the sphere is defined by

$G_{\varphi}(X,Y)=\mathrm{Re}(\xi,\eta).$ (3)

Perhaps surprisingly, the metric induced on the submanifold $M^{\sigma}_{3,3}$
of $CP^{L_{2}}$ with $2\sigma$ as a unit of length turns out to be the
ordinary Euclidean metric. In other words, the map
$\Omega:\mathbb{R}^{3}\times\mathbb{R}^{3}\longrightarrow CP^{L_{2}}$ is an
isometric embedding of the Euclidean space into the space of states. The
manifold $M^{\sigma}_{3,3}$ can be also furnished with a compatible linear
structure, making it isomorphic to the Euclidean space
$\mathbb{R}^{3}\times\mathbb{R}^{3}$.

A simple way to prove that Newtonian dynamics is a constrained Schrödinger
dynamics is by using the variational principle. The variation of the
functional

$S[\varphi]=\int\overline{\varphi}({\bf
x},t)\left[i\hbar\frac{\partial}{\partial t}-{\widehat{h}}\right]\varphi({\bf
x},t)d^{3}{\bf x}dt$ (4)

with the Hamiltonian ${\widehat{h}}=-\frac{\hbar^{2}}{2m}\Delta+V({\bf x},t)$
yields the Schrödinger equation for $\varphi$. For the states $\varphi$
constrained to the manifold $M^{\sigma}_{3,3}$, this functional reduces to the
classical action

$S=\int\left[{\bf p}\frac{d{\bf a}}{dt}-h({\bf p},{\bf a},t)\right]dt,$ (5)

where $h({\bf p},{\bf a},t)=\frac{{\bf p}^{2}}{2m}+V({\bf a},t)+const$ is the
Hamiltonian function for the particle. We used here the fact that for a
sufficiently small $\sigma$, the terms $V({\bf a},t)$ and $\int g^{2}_{{\bf
a},\sigma}({\bf x})V({\bf x},t)d^{3}{\bf x}$ are arbitrarily close to each
other. In fact, as $\sigma$ approaches $0$, the terms $g^{2}_{{\bf a},\sigma}$
form a delta sequence. It follows that the variation of the functional (4)
subject to the constraint that $\varphi$ belongs to $M^{\sigma}_{3,3}$ yields
Newtonian equations of motion.

The transition from Schrödinger to Newtonian dynamics can be expressed in
terms of the transition from quantum commutators to Poisson brackets. First of
all, a simple calculation with the above $h({\bf p},{\bf a},t)$ and
${\widehat{h}}$, and $\varphi$ in $M^{\sigma}_{3,3}$ yields the following:

$(\varphi,\tfrac{1}{i\hbar}[{\widehat{\bf x}},{\widehat{h}}]\varphi)=\\{{\bf
a},h\\},$ (6)

and

$(\varphi,\tfrac{1}{i\hbar}[{\widehat{\bf p}},{\widehat{h}}]\varphi)=\\{{\bf
p},h\\}.$ (7)

The brackets on the right side of (6) and (7) are the usual Poisson brackets.
Consider the linear vector fields $x_{\varphi}={\widehat{\bf x}}\varphi$ and
$p_{\varphi}={\widehat{\bf p}}\varphi$ on $L_{2}(\mathbb{R}^{3})$ associated
with the operators of position and momentum. For a given state $\varphi$, let
$\overline{\bf x}$ and $\overline{\bf p}$ be the expected values of the
operators on this state. The components $x_{\varphi\perp}=({\widehat{\bf
x}}-\overline{\bf x}I)\varphi$ and $p_{\varphi\perp}=({\widehat{\bf
p}}-\overline{\bf p}I)\varphi$ are tangent to the sphere $S^{L_{2}}$ at
$\varphi$ and orthogonal to the fibres $\\{\varphi\\}$ of the fibre bundle
$CP^{L_{2}}=S^{L_{2}}/S^{1}$. Therefore, they are tangent to the projective
space $CP^{L_{2}}$ itself. Moreover, the vector fields $x_{\varphi\perp}$ and
$p_{\varphi\perp}$ constrained to the submanifold $M^{\sigma}_{3,3}$ of
$CP^{L_{2}}$ are tangent to it. The integral curves of these fields are lines
of constant position and momentum that provide the coordinate grid in the
classical phase space $\mathbb{R}^{3}\times\mathbb{R}^{3}$. This can be seen
directly by comparing the right hand sides of the equations

$({\widehat{\bf x}}-\overline{\bf x}I)\varphi=({\bf x}-{\bf a})\varphi$ (8)

and

$({\widehat{\bf p}}-\overline{\bf p}I)\varphi=-\frac{i\hbar}{2\sigma^{2}}({\bf
x}-{\bf a})\varphi,$ (9)

valid for $\varphi$ in $M^{\sigma}_{3,3}$, with the gradients $\nabla_{\bf
a}\varphi$ and $\nabla_{\bf p}\varphi$. The latter gradients represent vectors
tangent to the lines of constant values of ${\bf p}$ and ${\bf a}$ through a
point of the projective manifold $M^{\sigma}_{3,3}$.

Using the Ehrenfest theorem written in terms of the vector fields
$x_{\varphi\perp}$ and $p_{\varphi\perp}$ and (6) and (7), we have for the
initial state $\varphi$ in $M^{\sigma}_{3,3}$ at $t=0$:

$2\mathrm{Re}\left(\frac{d\varphi}{dt},x_{\varphi\perp}\right)=\\{{\bf
a},h\\}=\frac{d{\bf a}}{dt}$ (10)

and

$2\mathrm{Re}\left(\frac{d\varphi}{dt},p_{\varphi\perp}\right)=\\{{\bf
p},h\\}=\frac{d{\bf p}}{dt}.$ (11)

So, at any point $({\bf a},{\bf p})$ of the phase space
$\mathbb{R}^{3}\times\mathbb{R}^{3}$, the derivatives $\frac{d{\bf a}}{dt}$
and $\frac{d{\bf p}}{dt}$ of a Newtonian motion are projections at the point
$\Omega({\bf a},{\bf p})$ of the velocity $\frac{d\varphi}{dt}$ of the
Schrödinger evolution onto the coordinate lines of the phase space submanifold
$M^{\sigma}_{3,3}$. This is a restatement of the variational principle for the
functional (4) with the state $\varphi$ constrained to manifold
$M^{\sigma}_{3,3}$.

Furthermore, as explained in [3], it is possible to decompose the initial
velocity of state of a particle in a potential $V$ with the initial state in
$M^{\sigma}_{3,3}$ into four physically meaningful orthogonal components. The
two components tangent to the phase space submanifold $M^{\sigma}_{3,3}$ were
already identified with the usual classical velocity ${\bf v}=d{\bf a}/dt$ and
acceleration ${\bf w}=d{\bf v}/dt=-\nabla V/m$ of the particle. The other two
components orthogonal to $M^{\sigma}_{3,3}$ are the phase velocity (projection
of $d\varphi/dt$ onto the unit vector $-i\varphi$) and the velocity of
spreading of the wave-packet (projection of $d\varphi/dt$ onto the unit vector
in the direction $id\varphi/d\sigma$). The norm of the total velocity
$d\varphi/dt$ at $t=0$, i.e., the speed of motion of the initial state
$\varphi$, is given by the formula

$\left\|\frac{d\varphi}{dt}\right\|^{2}=\frac{{\overline{E}}^{2}}{\hbar^{2}}+\frac{{\bf
v}^{2}}{4\sigma^{2}}+\frac{m^{2}{\bf
w}^{2}{\sigma}^{2}}{\hbar^{2}}+\frac{\hbar^{2}}{32\sigma^{4}m^{2}}.$ (12)

When the state is constrained to the manifold $M^{\sigma}_{3,3}$, the first
term (the phase velocity squared) and the last term (the velocity of spreading
squared) disappear. The only surviving terms are associated with the classical
velocity and acceleration and the motion follows the classical Newtonian
dynamics in the phase space.

Similar results are true for systems of $n$ particles. In this case, the
transition to classical Newtonian dynamics happens when the state of the
system is constrained to the phase space submanifold
$M^{\sigma}_{3n,3n}=M^{\sigma}_{3,3}\otimes...\otimes M^{\sigma}_{3,3}$ in the
space of states of the system. That is, the constrained state is the product
of states of the $n$ particles in the system with the state of each particle
constrained to a copy of $M^{\sigma}_{3,3}$. A more detailed analysis leads to
the following theorem [2]:

> (A) The Newtonian dynamics of an arbitrary mechanical system is the
> Schrödinger dynamics of that system with the state constrained to the
> classical phase space submanifold of the space of states of the system.
> Furthermore, the Schrödinger dynamics is the only unitary evolution that
> reduces under the constraint to the Newtonian one.

We note that the relationship between commutators and Poisson brackets
established in (6) and (7) generalizes directly to systems of $n$ particles.

The identification of the classical phase space of a particle with the
manifold $M^{\sigma}_{3,3}$ yields a useful metric relationship. Let
$M^{\sigma}_{3}$ denote the submanifold of $M^{\sigma}_{3,3}$ of the Gaussian
states $g_{{\bf
a},\sigma}=\left(\frac{1}{2\pi\sigma^{2}}\right)^{3/4}e^{-\tfrac{({\bf x}-{\bf
a})^{2}}{4\sigma^{2}}}$. The map $\omega:{\bf a}\longrightarrow g_{{\bf
a},\sigma}$, $\omega(\mathbb{R}^{3})=M^{\sigma}_{3}$ identifies the
submanifold $M^{\sigma}_{3}$ of $CP^{L_{2}}$ with the classical space
$\mathbb{R}^{3}$. If $\theta(g_{{\bf a},\sigma},g_{{\bf b},\sigma})$ denotes
the distance from $g_{{\bf a},\sigma}$ to $g_{{\bf b},\sigma}$ in
$M^{\sigma}_{3}$ in the Fubini-Study metric on the space of states
$CP^{L_{2}}$ and $({\bf a}-{\bf b})^{2}$ is the square of the Euclidean
distance between the corresponding points ${\bf a}$ and ${\bf b}$ in
$\mathbb{R}^{3}$, then

$e^{-\frac{({\bf a}-{\bf b})^{2}}{4\sigma^{2}}}=\cos^{2}\theta(g_{{\bf
a},\sigma},g_{{\bf b},\sigma}).$ (13)

Likewise, for arbitrary states $\varphi=g_{{\bf a},\sigma}e^{i{\bf p}{\bf
x}/\hbar}$ and $\psi=g_{{\bf b},\sigma}e^{i{\bf q}{\bf x}/\hbar}$ in the
classical phase space $M^{\sigma}_{3,3}$, the distances in the Euclidean phase
space and the space of states are related by

$e^{-\frac{({\bf a}-{\bf b})^{2}}{4\sigma^{2}}-\frac{({\bf p}-{\bf
q})^{2}}{\hbar^{2}/\sigma^{2}}}=\cos^{2}\theta(\varphi,\psi).$ (14)

The derivation of (13) and (14) is explained in the Appendix.

The relation (13) implies a connection between probability distributions of
multivariate random variables valued in the classical space and the space of
states. Namely, consider a random variable $\varphi$ with values in the space
of states $CP^{L_{2}}$. Because the classical space $M^{\sigma}_{3}$ is a
submanifold of $CP^{L_{2}}$, we can restrict $\varphi$ to take values in
$M^{\sigma}_{3}$. We then have the following result (see [2] and the
Appendix):

> (B) Suppose the conditional probability of $\varphi$ given that $\varphi$ is
> in the classical space is described by the normal distribution. Suppose the
> probability $P$ of transition between any two states depends only on the
> distance between the states. Then $P$ is given by the Born rule, i.e.,
> $P=|(\varphi,\varphi_{0})|^{2}$, where $\varphi_{0}$ is the initial and
> $\varphi$ is the observed state. The opposite is also true: For states in
> the classical space the Born rule yields the normal distribution of the
> observed position.

So, the normal probability distribution is the Born rule in disguise!

These results establish a deep connection between Newtonian and Schrödinger
dynamics and between the classical space and classical phase space, and the
submanifolds $M^{\sigma}_{3}$ and $M^{\sigma}_{3,3}$ of the space of states.
Based on these results we put forward the following embedding hypothesis:

> (EH) The constructed mathematical embedding of the classical space and
> classical phase space into the space of states and the resulting
> identification of Newtonian dynamics with the constrained Schrödinger
> dynamics are physical. That is, the classical space, phase space and the
> Newtonian dynamics of a system are not only fully derivable from, but also
> physically originate in the Schrödinger dynamics in the Hilbert space of
> states of the system.

The hypothesis asserts that the provided mathematical derivation of the
classical Newtonian dynamics from the Schrödinger dynamics of the system is
the true physical correspondence between classical and quantum systems.

To validate the hypothesis, we need to show that it is consistent with all
observed classical and quantum phenomena. This involves showing that all
observable classical phenomena can be derived from the corresponding quantum
phenomena by constraining the state to the classical phase space submanifold
in the space of states. This also involves physically explaining the
constraint itself. Now, an arbitrary deterministic classical motion is
described by the Newtonian equations of motion. The theorem ${\bf(A)}$
demonstrates that these equations follow from the Schrödinger equation with
the constraint. Therefore, given the constraint, the hypothesis is consistent
with an arbitrary deterministic classical motion. To validate the hypothesis,
it remains then to: (1) verify its consistency for the motion of particles
described in statistical mechanics, specifically, for the Brownian motion, and
(2) explain the origin of the constraint. The rest of the paper will be
dedicated to these two tasks.

Instead of diving into the subject of a quantum version of the Brownian
motion, note that a macroscopic Brownian particle in a medium can be
considered a classical chaotic system. For instance, the motion of the
particle through an appropriate lattice of round obstacles provides a
classical chaotic realization of the Brownian motion. A lively discussion of
the stochastic, deterministic chaotic, and regular characterizations of the
Brownian motion can be found in [12, 13, 14]. With the view that the Brownian
particle in a medium is a chaotic system comes the applicability of the BGS-
conjecture to the system [11]. It asserts that the Hamiltonian of the
corresponding quantum system can be represented by a random matrix.

Random matrices were originally introduced into quantum mechanics by Wigner
[10] in a study of excitation spectra of heavy nuclei. Wigner reasoned that
the motion of nucleons within the nucleus is so complex that the Hamiltonian
of the system can be modeled by a random matrix from an ensemble that respects
symmetries of the system but otherwise contains no additional information. It
was later discovered that correlations in the spectrum of random matrices
possess a remarkable universality in being applicable to a large number of
quantum systems. That includes nuclear, atomic, molecular and scattering
systems, chaotic systems with few degrees of freedom as well as complex
systems, such as solid state systems and lattice systems in the field theory.
This wealth of experimental evidence suggests that all quantum systems whose
classical counterpart is chaotic exhibit random matrix statistics. This is the
essence of the BGS-conjecture.

As mentioned, from the possibility to attribute chaotic character to Brownian
motion and the BGS conjecture it follows that the Hamiltonian of the quantum
analogue of the Brownian motion at any time is given by a random matrix. On
the physical grounds, we can also claim that the random matrices that
represent the Hamiltonian at two different moments of time must be
independent. This leads us to the following version of the BGS-conjecture:

> (RM) The quantum-mechanical analogue of the Brownian motion can be modeled
> by a random walk of state in the space of states of the system. The steps of
> the random walk without drift satisfy the Schrödinger equation with the
> Hamiltonian represented by a random matrix from the Gaussian unitary
> ensemble (GUE). The matrices representing the Hamiltonian at two different
> times are independent and belong to the same ensemble.

Note that the Schrödinger equation with the Hamiltonian represented by a
random matrix describes evolution of the state and not of the density matrix
for the system. In that sense, it is analogous to the Langevin equation for
the position of a Brownian particle rather than the diffusion equation for the
probability density function.

First of all, let us prove that the random walk in ${\bf(RM)}$ yields the
Brownian motion in $\mathbb{R}^{3}$. More precisely,

> (C) The random walk described in ${\bf(RM)}$ but conditioned to stay on the
> submanifold $M^{\sigma}_{3}$ in the space of states yields a random walk on
> $\mathbb{R}^{3}$ that approximates the Brownian motion of a particle in a
> medium.

In fact, a general Schrödinger evolution with Hamiltonian ${\widehat{h}}$ can
be thought of as a sequence of steps connecting the points
$\varphi_{t_{0}},\varphi_{t_{1}},...$ in the space of states. For small time
intervals $\Delta t=t_{k}-t_{k-1}$, the state $\varphi_{t_{N}}$ at time
$t_{N}$ is given by the time ordered product

$\varphi_{t_{N}}=e^{-\frac{i}{\hbar}{\widehat{h}}(t_{N})\Delta
t}e^{-\frac{i}{\hbar}{\widehat{h}}(t_{N-1})\Delta
t}...e^{-\frac{i}{\hbar}{\widehat{h}}(t_{1})\Delta t}\varphi_{t_{0}}.$ (15)

Suppose the evolution of the state $\varphi$ of a particle is constrained to
the classical space submanifold $M^{\sigma}_{3}$. The points
$\varphi_{t_{0}},\varphi_{t_{1}},...$ belong then to the submanifold
$M^{\sigma}_{3}$ and the steps can be identified with translations in the
classical space. This is to say that for each $k$, the operator
${\widehat{h}}(t_{k})$ acts as the generator of translation by a vector
${\bf\xi}_{k}$ in $\mathbb{R}^{3}$, so that
${\widehat{h}}(t_{k})={\bf\xi}_{k}{\widehat{\bf p}}$, where ${\widehat{\bf
p}}$ is the momentum operator. Because all operators of translation commute
with each other, the equation (15) yields the following expression:

$\varphi_{t_{N}}({\bf x})=\varphi_{t_{0}}(x-{\bf\xi}_{1}\Delta
t-{\bf\xi}_{2}\Delta t-...-{\bf\xi}_{N}\Delta t).$ (16)

That is, the initial state is simply translated by the vector

${\bf d}=\sum^{N}_{k=1}{\bf\xi}_{k}\Delta t$ (17)

in $\mathbb{R}^{3}$. Now, the probability distribution of steps
$-\frac{i}{\hbar}{\widehat{h}}(t_{k+1})\varphi_{t_{k}}$ in the tangent space
$T_{\varphi_{k}}M^{\sigma}_{3}$ must be the conditional probability
distribution of steps for the Hamiltonian satisfying ${\bf(RM)}$ under the
condition that the steps take place in $T_{\varphi_{k}}M^{\sigma}_{3}$. From
the properties of the random matrix, it follows that ${\bf\xi}_{k}$ are
independent and identically normally distributed random vectors, so that the
equation (17) defines a random walk with Gaussian steps in $\mathbb{R}^{3}$.
This is known to approximate the Brownian motion in $\mathbb{R}^{3}$ and yield
the normal distribution of the position vector ${\bf d}$ at any time $t$,
proving the claim.

$\hfill\square$

By the extension of a random walk on $\mathbb{R}^{3}$ to the space of states
$CP^{L_{2}}$, we understand a walk in $CP^{L_{2}}$ that satisfies ${\bf(RM)}$
and that reduces to the original random walk on $\mathbb{R}^{3}$ when
conditioned to stay on $M^{\sigma}_{3}$. From ${\bf(C)}$, we know that an
extension of the walk with Gaussian steps on $\mathbb{R}^{3}$ exists. We claim
that such an extension is unique. In fact, because the random walk conditioned
to stay on $M^{\sigma}_{3}$ must be Gaussian, the entries of the matrix of the
Hamiltonian associated with the directions tangent to $M^{\sigma}_{3}$ at a
point are distributed normally. But the probability distribution of a single
entry of the matrix defines the corresponding Gaussian unitary ensemble.
Because of that, the random walk with Gaussian steps in $\mathbb{R}^{3}$
defines a unique random walk satisfying ${\bf(RM)}$. In what follows, the
random walk of the state of a particle in ${\bf(RM)}$ will always be assumed
to be this unique extension of the random walk with Gaussian steps that
approximate the Brownian motion of the particle. We are free to make this
choice as long as the resulting extension satisfies our needs. At the same
time, this choice makes sense physically, because a displacement of the
particle in $\mathbb{R}^{3}$ results in the like displacement of its state in
$M^{\sigma}_{3}$, tying the two motions together.

The theorems ${\bf(A)}$ and ${\bf(C)}$ prove the consistency of the hypothesis
${\bf(EH)}$ for the deterministic Newtonian dynamics and the Brownian motion.
Our goal now is to explain the origin of the constraint to the classical space
and phase space submanifolds in the Hilbert space of states. For this, we need
to establish properties of the random walk on the space of states in the
conjecture ${\bf(RM)}$ and establish its applicability to the process of
measurement.

> (D) The probability distribution of steps of the random walk specified in
> ${\bf(RM)}$ is isotropic and homogeneous. That is, for all initial states
> $\\{\varphi\\}$ in the space of states $CP^{H}$, the vector
> $d\varphi=-\frac{i}{\hbar}{\widehat{h}}\varphi dt$ is a normal random vector
> in the tangent space $T_{\\{\varphi\\}}CP^{H}$, with spherical distribution.
> The probability of reaching a certain state during the walk depends only on
> the distance between the initial and the final state.

To prove, note that because for any $t$ the matrix of ${\widehat{h}}$ is in
GUE, the probability density function $P({\widehat{h}})$ of ${\widehat{h}}$ is
invariant with respect to conjugations by unitary transformations. That is,
$P(U^{-1}{\widehat{h}}U)=P({\widehat{h}})$ for all unitary transformations $U$
acting in the Hilbert space of states. Also, for all unitary transformations
$U$ that leave $\\{\varphi\\}$ unchanged and therefore all $U$ that act in the
tangent space $T_{\\{\varphi\\}}CP^{H}$, we have

$(U^{-1}{\widehat{h}}U\varphi,v)=({\widehat{h}}U\varphi,Uv)=({\widehat{h}}\varphi,Uv),$
(18)

where $v$ is a unit vector in $T_{\\{\varphi\\}}CP^{H}$. It follows that

$P({\widehat{h}}\varphi,v)=P({\widehat{h}}\varphi,Uv),$ (19)

where $P$ is the probability density of the components of
${\widehat{h}}\varphi$ in the given directions. By a proper choice of $U$, we
can make $Uv$ to be an arbitrary unit vector in $T_{\\{\varphi\\}}CP^{H}$,
proving the isotropy of the distribution. On the other hand, for all unitary
operators $V$ in $H$ and a unit vector $v$ in $T_{\\{\varphi\\}}CP^{H}$, we
have

$P({\widehat{h}}\varphi,v)=P(V^{-1}{\widehat{h}}V\varphi,v)=P({\widehat{h}}V\varphi,Vv).$
(20)

Because $V\varphi$ is an arbitrary state and $Vv$ is in the tangent space
$T_{\\{V\varphi\\}}CP^{H}$, we conclude with the help of (18) that the
probability density function is independent of the initial state of the
system, proving the homogeneity of the distribution. The components of the
vector ${\widehat{h}}\varphi dt$ are independent, because the entries of the
matrix of ${\widehat{h}}$ are independent. It follows that
$-\frac{i}{\hbar}{\widehat{h}}\varphi dt$ is a normal random vector with
spherical distribution. Finally, because different steps of the walk are
independent, the probability of reaching a certain state during the walk may
depend only on the distance between the initial and the final state.

$\hfill\square$

Let us show now that the version ${\bf{(RM)}}$ of the BGS-conjecture provides
a consistent approach to the process of observation of the position of a
single particle in the classical and quantum mechanics alike.

> (E) The Schrödinger evolution with the Hamiltonian that satisfies
> ${\bf(RM)}$ models the measurement of the position of a particle in
> classical and quantum physics. Namely, under the constraint that the random
> walk that approximates the evolution of state takes place on the submanifold
> $M^{\sigma}_{3}$, the probability distribution of the position random vector
> is normal. Without the constraint, the probability to find the initial state
> $\varphi_{0}$ at a point $\varphi$ is given by the Born rule. The transition
> from the initial to final state is time-irreversible.

In fact, according to ${\bf(C)}$, the random walk that satisfies ${\bf(RM)}$
but is conditioned to stay on $M^{\sigma}_{3}$ is the random walk with
Gaussian steps on $\mathbb{R}^{3}$. The latter random walk considered over the
time interval of observation yields the normal distribution of the position
random vector. Note that the normal distribution of the position agrees with
observations in the macro world. It is also consistent with the central limit
theorem applied to describe the cumulative effect of uncontrollable
fluctuations from the mean in a series of measurement outcomes. On the other
hand, according to ${\bf(D)}$, the probability $P$ of reaching a state
$\varphi$ from the initial state $\varphi_{0}$ by means of a unconstrained
random walk that satisfies ${\bf(RM)}$ depends only on the distance between
the states. From ${\bf(B)}$, it then follows that $P$ is given by the Born
rule: $P=|(\varphi,\varphi_{0})|^{2}$. Also, the choice of the Gaussian
unitary ensemble corresponds to time-irreversible systems [15]. It follows
that the model satisfies the basic properties of observation of the position
in classical and quantum mechanics.

$\hfill\square$

Quantum mechanics does not explain how and when the deterministic Schrödinger
evolution is replaced by probabilistic evolution, whose outcomes obey the Born
rule. It does not explain why macroscopic bodies are never observed in
superpositions of position eigenstates. The wave-function collapse models
[19],[20] aim to answer these questions by introducing a non-linear stochastic
modification of the Schrödinger equation. It is assumed that the modified
equation must be non-linear to be able to suppress superpositions of states.
It is also assumed that the modified equation must vary with a change in the
observed quantity. However, here the Born rule was derived from the linear
Schrödinger equation with the Hamiltonian satisfying ${\bf(RM)}$. Moreover, it
was derived for all observed quantities at once, without needed to change the
equation of motion. Is there a contradiction?

The modified Schrödinger equations in the collapse models make the state of
the system converge to an eigenstate of the measured quantity, usually the
position of a particle. An equation like that must break an arbitrary initial
superposition of states, so it must be non-linear. Under the Schrödinger
evolution with the Hamiltonian satisfying the conjecture ${\bf(RM)}$, the
state of the system does not converge to a position eigenstate. The equation
does not suppress superpositions. It makes the state wander around the space
of states and ensures that the probability of reaching a particular
neighborhood in the space is given by the Born rule. That explains why the two
approaches do not contradict each other.

We are now ready to investigate why the state of a macroscopic body is
confined to the classical space submanifold of the space of states. The key to
this is that the Brownian motion and the motion of state under a measurement
are now derived from the same dynamics. We know that large particles do not
experience a Brownian motion. That is because the total force acting on any
such particle from the particles of the medium is nearly zero. As a result, in
the absence of an external potential, the particle remains at rest in the
medium. A similar mechanism explains why the state of a particle may be
confined to the submanifold $M^{\sigma}_{3}$:

> (F) Suppose the Hamiltonian of a particle in the natural surroundings
> satisfies the conjecture ${\bf(RM)}$. Then, the state of the particle is
> constrained to the submanifold $M^{\sigma}_{3}$ precisely when the induced
> Brownian motion of the particle described in ${\bf(C)}$ vanishes. More
> precisely, under ${\bf(RM)}$, the boundary between the quantum and the
> classical occurs for particles that satisfy the following two conditions.
> The particles must be sufficiently large in size so that their Brownian
> motion in the natural surroundings is observable. At the same time, the
> particles must not be too large, when the Brownian motion trivializes, the
> state initially in $M^{\sigma}_{3}$ becomes constrained to $M^{\sigma}_{3}$
> and the Schrödinger evolution becomes Newtonian motion.

Recall that the spaces $M^{\sigma}_{3}$ and $\mathbb{R}^{3}$ are metrically
identical. Recall also that the random walk in ${\bf(RM)}$ is the unique
extension of the random walk in $\mathbb{R}^{3}$ that approximates the
Brownian motion of the particle in the surroundings. Consider a particle whose
initial state is in $M^{\sigma}_{3}$. Suppose the particle is sufficiently
large in size so that the Brownian motion of the particle in the natural
surroundings is negligible. Therefore, the motion of state in the directions
tangent to $M^{\sigma}_{3}$ is negligible. From ${\bf(D)}$, we know that for
the Hamiltonian satisfying ${\bf(RM)}$, the probability distribution of steps
of the random walk in the space of states of the particle is isotropic.
Because the probability of steps in the directions tangent to $M^{\sigma}_{3}$
vanishes, the same holds true for any other direction tangent to the space of
states at the same point. Therefore, in the absence of an external potential,
the motion of state of the particle in the space of states is trivial. When an
external potential is applied to the particle, the two middle terms in the
decomposition (12) may appear. However, these terms can only contribute to the
motion of state in the direction tangent to the classical phase space
submanifold $M^{\sigma}_{3,3}$. Therefore, the state will remain constrained
to the submanifold. The theorem ${\bf{(A)}}$ asserts then that the particle in
such a state will move in accord with the Newtonian dynamics.

On the contrary, suppose the particle is sufficiently small, but not too
small, so that the interaction between the particle and the surroundings
cannot be ignored and results in a noticeable Brownian motion of the particle.
By the isotropy of the probability distribution of steps of the random walk of
state, a displacement of the particle away from $M^{\sigma}_{3}$ is then
equally likely. Such a displacement would mean that the particle is now in a
superposition of states of different positions. If ${\bf a}$ is the initial
position and ${\bf l}$ is the observed displacement of the particle in
$M^{\sigma}_{3}$ during the measurement, then the states $g_{{\bf a},\sigma}$
and $g_{{\bf a}+{\bf l},\sigma}$ are distinguishable in the experiment. It
means that the superpositions of these states can be observed as well,
indicating that we are dealing with a quantum system.

$\hfill\square$

We see how the properties of the Hamiltonian in the conjecture ${\bf(RM)}$ are
responsible for the fact that the state of a macroscopic particle driven by
the Schrödinger dynamics with this Hamiltonian is constrained to the classical
phase space manifold. One might be tempted to say that the resulting
“freezing” of the state is in agreement with the quantum Zeno effect for the
particle whose position is continuously measured. However, the essential
difference is that the result is derived here from the unitary evolution
generated by a random Hamiltonian without ever needing to involve the
projection operators.

As an example, in the Appendix we find the displacement of a particle of
radius $1mm$ in still air in normal conditions. Assuming the particle’s
position is measured by scattering and observing visible light, the
displacement during the time interval of the measurement is found to be of the
order of $10^{-12}m$. The corresponding displacement of state in the space of
states is of the order of $10^{-7}rad$. These displacements are too small to
be observed in the described experiment. The state of the particle remains
“frozen”. The particle is constrained to $M^{\sigma}_{3}$ and no
superpositions of states of a different position of the particle can be
observed.

As already discussed, the phase space of $n$ particles in the classical space
$\mathbb{R}^{3}$ can be identified with the $6n$-dimensional submanifold
$M^{\sigma}_{3n,3n}=M^{\sigma}_{3,3}\otimes M^{\sigma}_{3,3}\otimes...\otimes
M^{\sigma}_{3,3}$ in the space of states of the system. By ${\bf(A)}$, the
Schrödinger dynamics for the system with the state constrained to
$M^{\sigma}_{3n,3n}$ is equivalent to its Newtonian dynamics. We know that
sufficiently small macroscopic particles in a medium exhibit a Brownian
motion. To establish the underlying quantum process, we need to find the
Hilbert space of states of the system and to apply the conjecture ${\bf(RM)}$.
For instance, whenever the interaction between the particles can be neglected,
the Brownian motion of each particle can be considered independently. In this
case ${\bf(RM)}$ is applicable to each particle and the previously obtained
results apply. In the opposite case, when the particles interact strongly and
form a rigid body, the body as a whole undergoes a Brownian motion in the
medium.

As an example of the latter case, consider two particles connected by a
weightless rigid rod of length $\Delta$. Then the state of the pair in
$M^{\sigma}_{3}\otimes M^{\sigma}_{3}$ is given by
${\widetilde{\delta}}^{3}_{\bf a}\otimes{\widetilde{\delta}}^{3}_{{\bf
a}+{\bf\Delta}}$, where ${\bf\Delta}$ is a vector of length $\Delta$ in
$\mathbb{R}^{3}$. The Hilbert space $H$ of states for the system is a
completion of the space of all linear combinations of such states. Under
${\bf(RM)}$, the state of the measured system undergoes a random walk in the
space $CP^{H}$ of the obtained entangled states. The positions of both
particles are known when the position of one of them is known, as the path of
the walk passes in this case through a point of the manifold
$M^{\sigma}_{3}\otimes M^{\sigma}_{3}$. As before, the probability of passing
through a particular point is given by the Born rule.

Consider now a system consisting of a microscopic particle $P$ and a
macroscopic device $D$ capable of measuring the position of $P$. The particle
and the device form a two-particle system whose state belongs to the product
Hilbert space $H=H_{P}\otimes H_{D}$. Consider the motion of the corresponding
classical system. Let say, $D$ is a cloud chamber, small enough to be
considered a material point, when treated classically. The macroscopic device
interacts with the surroundings. It also interacts with the particle. However,
in classical physics, the effect of the particle on the device can be
neglected. The interaction of the device with the surroundings results in a
Brownian motion of the device. When the device is sufficiently large, its
Brownian motion is trivial and the device is at rest in the lab system. By
applying ${\bf(RM)}$ to the device, we conclude that the state of the device
positioned initially in $M^{\sigma}_{3}$ can be treated independently of the
state of the particle and is at rest in the space of state $CP^{H_{D}}$. On
the other hand, a small macroscopic particle placed in the medium of the
chamber would undergo a Brownian motion. By applying ${\bf(RM)}$ to this case,
we conclude that the state of the corresponding quantum system will perform a
random walk in the space $CP^{H_{P}}$. According to ${\bf(E)}$, the
probability to find the particle at a particular point of $\mathbb{R}^{3}$
during the walk is given by the Born rule.

By the previous discussion, the deterministic Schrödinger evolution and the
random walk that satisfies ${\bf(RM)}$ are unique extensions of their
classical counterparts from a classical space submanifold to the space of
states of the system. We conclude that the particle-device system initially in
a product state must remain in a product state during the evolution.
Furthermore, according to ${\bf(F)}$, the state of the device initially
positioned in $M^{\sigma}_{3}$ will be confined to this submanifold. In
particular, the theorem ${\bf(A)}$ applies and ensures that in the presence of
an external potential, the device will be evolving in accord with the
Newtonian dynamics. Under these conditions, the system remains in a product
state with both factors being able to change. In particular, the position of
the particle can be mechanically recorded with no entanglement between the
particle and the device ever appearing or needing to appear in the process. To
summarize:

> (G) Suppose the initial state of the system consisting of a microscopic and
> a macroscopic particle in the natural surroundings is a product state with
> the state of the macroscopic particle in the manifold $M^{\sigma}_{3}$.
> Suppose also that the conjecture ${\bf(RM)}$ holds true. Then, during the
> evolution, the state of the system remains in a product form with the state
> of macroscopic particle confined to $M^{\sigma}_{3}$ and evolving
> classically.

Theorem ${\bf(G)}$ provides an explanation of the issue of an apparent
inconsistency of the description of quantum mechanics by macroscopic
observers. The issue was pointed out first by Wigner [16] and has picked up a
significant interest in recent times. The reason for the contradiction in the
accounts of Wigner and the friends in the Wigner’s friend-type thought
experiments is the assumption that a macroscopic observer or the whole lab may
exist in a superposition of classical states. Provided the conjecture
${\bf(RM)}$ for systems in the natural surroundings holds true, the evolution
of the state of the observers and the lab in the space of states can only
happen within the classical space submanifold. The state of the measured
system and the observer or the lab remains the product state. This invalidates
the no-go theorem in [17] and related results claiming the inconsistency of
accounts by the observers. Note that the experiment reported in [18], in which
the photons played the role of an observer does not change this conclusion.
The reality as “described” by microscopic particles is different than the one
the macroscopic observers are aware of.

Note also that the random matrix description of the interaction between
macroscopic bodies and the surroundings allowed us to bypass the decoherence-
based approach to the transition from quantum to classical. The modifications
of the Schrödinger equation with the goal of accounting for the measurement
results and explaining the transition to classicality [19, 20] become
unnecessary as well. The single outcome of a measurement, the probability of
the outcome and the transition from quantum to classical all have now a simple
explanation, rooted in the Schrödinger equation, the embedding of the
classical space into the space of states and the conjecture ${\bf(RM)}$. Let
us also remark that the elegance and the scope of the proposed solution to the
issue of the relationship between the quantum and classical systems suggests
that natural phenomena are happening in the space of states and not in the
classical space. If this is the case, the space of states becomes a new arena
for all physical processes with the classical space of macroscopic objects as
a submanifold of the space of states.

##  2 Appendix

The isomorphism of the classical space and classical phase space and the
submanifolds $M^{\sigma}_{3}$ and $M^{\sigma}_{3,3}$ in the space of states
$CP^{L_{2}}$ is explained in [3, 4]. The decomposition (12) is obtained by
projecting the velocity $\frac{d\varphi}{dt}$ of state under the Schrödinger
evolution with an arbitrary Hamiltonian
${\widehat{h}}=-\frac{\hbar^{2}}{2m}\Delta+V({\bf x},t)$ at a point
$\varphi=g_{{\bf a},\sigma}e^{i{\bf p}{\bf x}/\hbar}$ in the classical phase
space $M^{\sigma}_{3,3}$ onto an orthonormal set of vectors specified by
changing the parameters, ${\bf a},{\bf p},\sigma$ that define $\\{\varphi\\}$
and the phase parameter $\theta$ of a possible constant phase factor
$e^{-i\theta}$ of $\varphi$. Calculation of the classical space components of
$\frac{d\varphi}{dt}$ at an arbitrary point $\varphi$ in $M^{\sigma}_{3,3}$
yields

$\left.\mathrm{Re}\left(\frac{d\varphi}{dt},-\widehat{\frac{\partial\varphi}{\partial
a^{\alpha}}}\right)\right|_{t=0}=\frac{p^{\alpha}}{2m\sigma},$ (21)

where the hat here and in other calculations denotes normalization. For the
momentum space components of $\frac{d\varphi}{dt}$ at $\varphi$ we similarly
obtain, assuming that $\sigma$ is small enough to make the linear
approximation of $V({\bf x})$ valid:

$\left.\mathrm{Re}\left(\frac{d\varphi}{dt},\widehat{\frac{\partial\varphi}{\partial
p^{\alpha}}}\right)\right|_{t=0}=\frac{mw^{\alpha}\sigma}{\hbar},\quad{\text{w}here}\quad
mw^{\alpha}=-\left.\frac{\partial V({\bf x})}{\partial
x^{\alpha}}\right|_{{\bf x}={\bf a}}.$ (22)

The components (21) and (22) are tangent to $M^{\sigma}_{3,3}$ and orthogonal
to the fibre $\\{\varphi\\}$. The component of the velocity
$\frac{d\varphi}{dt}$ due to change in $\sigma$ (spreading) is orthogonal to
the phase space $M^{\sigma}_{3,3}$ and the fibre $\\{\varphi\\}$, and is equal
to

$\mathrm{Re}\left(\frac{d\varphi}{dt},i\widehat{\frac{d\varphi}{d\sigma}}\right)=\frac{\sqrt{2}\hbar}{8\sigma^{2}m}.$
(23)

The component of the velocity parallel to the fibre $\\{\varphi\\}$ is the
expected value of energy divided by $\hbar$:

$\mathrm{Re}\left(\frac{d\varphi}{dt},-\widehat{\frac{d\varphi}{d\theta}}\right)=\frac{1}{\hbar}(i{\widehat{h}}\varphi,i\varphi)=\frac{{\overline{E}}}{\hbar}.$
(24)

Calculation of the norm of
$\frac{d\varphi}{dt}=\frac{i}{\hbar}{\widehat{h}}\varphi$ at $t=0$ gives

$\left\|\frac{d\varphi}{dt}\right\|^{2}=\frac{{\overline{E}}^{2}}{\hbar^{2}}+\frac{{\bf
p}^{2}}{4m^{2}\sigma^{2}}+\frac{m^{2}{\bf
w}^{2}{\sigma}^{2}}{\hbar^{2}}+\frac{\hbar^{2}}{32\sigma^{4}m^{2}},$ (25)

which is the sum of squares of the found components. This completes a
decomposition of the velocity of state at any point $\varphi$ in
$M^{\sigma}_{3,3}$.

The metric relationship (13) follows from the inner product of two states in
$M^{\sigma}_{3}$:

$(g_{{\bf a},\sigma},g_{{\bf b},\sigma})=e^{-\frac{({\bf a}-{\bf
b})^{2}}{8\sigma^{2}}}.$ (26)

This expression squared is equal to the right hand side of (13) by the
definition of the Fubini-Study distance between states in $CP^{L_{2}}$. The
result (14) is obtained in a similar way by evaluating the Fourier transform
of a Gaussian function along the way.

The estimate of the diffusion coefficient of a macroscopic spherical particle
of radius $1mm$ in still air in normal conditions is based on the Stokes-
Einstein relationship

$D=\frac{k_{B}T}{6\pi\eta r},$ (27)

where $D$ is the diffusion coefficient, $r$ is the radius of the particle,
$\eta$ is the dynamic viscosity, $T$ is temperature of the medium and $k_{B}$
is the Boltzmann constant. Using the room temperature and the known value of
dynamic viscosity $\eta\sim 10^{-5}N\cdot s/m^{2}$, we get $D\sim
10^{-12}m^{2}/s$. The variance of the $x$-coordinate of position of the
particle is given by $\overline{x^{2}}=2Dt$. If we scatter visible light off
the particle to determine its position, the time interval of observation can
be estimated to be as short as $10^{-13}s$. This amounts to the displacement
of the order of $10^{-12}m$. The accuracy of measurement is limited by the
wavelength $\lambda\sim 10^{-5}m$. The Fubini-Study distance between Gaussian
states that are $10^{-12}m$ apart in $M^{\sigma}_{3}$ with $\sigma\sim
10^{-5}m$ is calculated via (13) and is about $10^{-7}rad$.

## References

## References

  * [1] Kryukov, A. On the motion of macroscopic bodies in quantum theory. J. Math. Phys. 58, 082103 (2017).
  * [2] Kryukov, A. On observation of position in quantum theory. J. Math. Phys. 59 052103 (2018).
  * [3] Kryukov, A. Mathematics of the classical and the quantum. J. Math. Phys. 61 082101 (2020).
  * [4] Kryukov, A. Linear algebra and differential geometry on abstract Hilbert space. Int. J. Math. Math. Sci. 2005 2241 (2005).
  * [5] Kryukov, A. Geometric derivation of quantum uncertainty. Phys. Lett. A 370 419-422 (2007).
  * [6] Ullersma, P. An exactly solvable model for Brownian motion. Physica 32 74-96 (1966).
  * [7] Zwanzig, R. Nonequilibrium Statistical Mechanics. (Oxford University Press, New York, 2001).
  * [8] Gaspard, P. Chaos, scattering and Statistical Mechanics. (Cambridge University Press, Cambridge, 1998).
  * [9] Szasz, D., editor. Hard Ball Systems and the Lorentz Gas. (Springer-Verlag, Berlin, 2000).
  * [10] Wigner, E. P. On the statistical distribution of the widths and spacings of nuclear resonance levels. Proc. Cambridge Philos. Soc. 47 790 (1951).
  * [11] Bohigas, O., Giannoni, M. J. and Schmit, C. Characterization of Chaotic Quantum Spectra and Universality of Level Fluctuation Laws. Physical Review Letters 52 1 (1984).
  * [12] Gaspard, P., Briggs, M. E., Francis, M. K., Sengers, J. V., Gammon, R. W., Dorfman, J. R. and Calabrese, R. V. Experimental evidence for microscopic chaos. Nature 394 865-868 (1998).
  * [13] Grassberger, P. and Schreiber, T. Microscopic chaos from brownian motion? Nature 401 875-876 (1999).
  * [14] Cecconi, F., Cencini, M., Falcioni, M. and Vulpiani, A. Brownian motion and diffusion: from stochastic processes to chaos and beyond. Chaos 15(2) 26102 (2005).
  * [15] Wigner, E.P. Group theory and its application to quantum mechanics of atomic spectra. (Academic Press, New York, 1959).
  * [16] Wigner, E.P. Remarks on the mind-body question. In: Good,I.J., editor. The Scientist Speculates, (Heinemann, London, 1961)
  * [17] Frauchiger, D. and Renner, R. Quantum theory cannot consistently describe the use of itself. Nature Communications 9 3711 (2018).
  * [18] Proietti, M., Pickston, A., Graffitti, F., Barrow, P., Kundys, D., Branciard, C., Ringbauer, M. and Fedrizzi, A. Experimental test of local observer independence. Science Advances 5 9832 (2019).
  * [19] Bassi, A. and Ghirardi, G. Dynamical Reduction Models. Physics Reports 379 257 (2003).
  * [20] Bassi, A., Lochan, K., Satin, S., Singh, T. and Ulbricht, H. Models of wave-function collapse, underlying theories, and experimental tests. Rev. Mod. Phys. 85 471 (2013).