A violent error via quil_to_native_quil on ARM64 #842
Comments
|
Would you be able to share a minimal program that reproduces the issue? |
|
@notmgsk It's not exactly minimal, but here's the circuit that caused the above error in a Edit: though of course now that I've posted this, it seems to compile without issue 🤔 |
|
Ok, I've got another failing example, this time on Aspen-11. Some traceback is snipped for brevity &/or IP concerns: Python traceback:
Offending (again, long, sorry, I'll try to shorten this) program: |
|
It's not happening reliably, but the following (deleted as much as I could from the above that didn't touch qubits 32, 45, or 46) crashes most of the time with the same Here's the |
|
Is this happening on ARM only? I find that this function is pretty sensitive to the underlying numerical libraries. I myself was experimenting with changing the eigenvalue function and was getting the same failures when doing so. Unfortunately I didn't dig deep enough to truly find out why. One thing we could try to do is replace this function wholesale with a more direct approach a la #841. |
|
So far it is only failing on ARM; the compilation succeeds on x86. I did see #841 on the issues list when I started this post and thought I'd been scooped! Regarding
Would that be a long and laborious process? |
I think the labor would be measured in hours. It would consist of:
|
This was referenced Sep 6, 2022
|
For another piece of information, here's what happened when trying to run |
In summary, I don't know what's going onI'm well and truly flummoxed at this point. Per a suggestion out-of-band, I attempted running the diagonalization routine and comparing both the # attempts and the residuals (largest absolute values of the off-diagonal entries of the matrix
|
| Platform | LAPACK | Installed via |
|---|---|---|
| Apple Laptop (Intel) | liblapack.3.10.1.dylib |
brew install lapack |
| x86 | liblapack-dev version 3.7.1-4ubuntu1 |
apt-get install liblapack-dev |
| ARM | liblapack-dev version 3.9.0-1build1 |
apt-get install liblapack-dev |
Python attempts:
Python residuals:
Lisp attempts:
Lisp residuals:
Supporting files
Python script gather_data.py
from typing import NamedTuple
import numpy as np
import pandas as pd
from scipy.stats import unitary_group
import typer
class _Decomp(NamedTuple):
attempts: int
eigvals: np.ndarray
eigvecs: np.ndarray
def _orthogonal_diagonalization(m: np.ndarray, rng: np.random.Generator) -> _Decomp:
"""github.com/gecrooks/quantumflow-dev/blob/master/quantumflow/decompositions.py#L325"""
for i in range(16):
c = rng.uniform(0, 1)
matrix = c * m.real + (1 - c) * m.imag
_, eigvecs = np.linalg.eigh(matrix)
eigvals = eigvecs.T @ m @ eigvecs
if np.allclose(m, eigvecs @ eigvals @ eigvecs.T):
return _Decomp(attempts=i, eigvals=eigvals, eigvecs=eigvecs)
raise np.linalg.LinAlgError("Matrix is not orthogonally diagonalizable")
class _Datum(NamedTuple):
attempts: int
max_off_diagonal: float
max_difference: float
def _collect(u: np.ndarray, rng: np.random.Generator) -> _Datum:
uut = u @ u.T
decomp = _orthogonal_diagonalization(uut, rng)
d, o = decomp.eigvals, decomp.eigvecs
m = o.T @ uut @ o
return _Datum(
attempts=decomp.attempts,
max_off_diagonal=np.abs(m - np.diag(np.diag(m))).max(),
max_difference=np.abs(o @ d @ o.T - uut).max(),
)
def collect_specifics(_iterations: int, rng: np.random.Generator) -> pd.Series:
residuals = []
for unitary in [
np.array([[1, 0, 0, 0], [0, 0, 0, 1], [0, 0, 1, 0], [0, 1, 0, 0]]),
np.diag([1, 1, 1, 1j]),
np.kron(np.eye(2), 1 / np.sqrt(2) * np.array([[1, -1j], [-1j, 1]])),
np.array(
[
[-0.965 + 0.133j, 0.093 + 0.000j, -0.000 + 0.000j, -0.000 + 0.206j],
[0.093 + 0.000j, 0.965 + 0.133j, 0.000 - 0.206j, -0.000 + 0.000j],
[-0.000 + 0.000j, 0.000 - 0.206j, 0.965 - 0.133j, -0.093 + 0.000j],
[0.000 + 0.206j, -0.000 + 0.000j, -0.093 + 0.000j, -0.965 - 0.133j],
]
),
]:
residuals.append(_collect(unitary, rng))
return pd.Series(residuals)
def collect_su4(iterations: int, rng: np.random.Generator) -> pd.Series:
residuals = []
gen = unitary_group(dim=4, seed=rng.integers(1 << 32))
for _ in range(iterations):
residuals.append(_collect(gen.rvs(), rng))
return pd.Series(residuals)
def collect_su2_x_su2(iterations: int, rng: np.random.Generator) -> pd.Series:
residuals = []
gen = unitary_group(dim=2, seed=rng.integers(1 << 32))
for _ in range(iterations):
residuals.append(_collect(np.kron(gen.rvs(), gen.rvs()), rng))
return pd.Series(residuals)
def main(
iterations: int = typer.Option(10_000, help="Number of rounds to run"),
seed: int = typer.Option(96692877, help="Seed for PRNG (random.org)"),
platform: str = typer.Option("laptop", help="Platform collected on"),
):
rng = np.random.default_rng(seed)
collections = []
for (name, func) in [
("Specifics", collect_specifics),
("SU(4)", collect_su4),
("SU(2)xSU(2)", collect_su2_x_su2),
]:
for d in func(iterations, rng):
collections.append(
{
"name": name,
"attempts": d.attempts,
"max_off_diagonal": d.max_off_diagonal,
"max_difference": d.max_difference,
}
)
df = pd.DataFrame(collections)
df["platform"] = platform
df.to_feather(f"{platform}_python.arrow")
if __name__ == "__main__":
typer.run(main)Lisp script gather-details.lisp
;;;; sbcl --load gather-data.lisp --eval "(sb-ext:save-lisp-and-die \"gather-data\" :executable t :save-runtime-options t :toplevel 'gather-data:toplevel)"
;;;; ./gather-data -p PLATFORM
;;;;
;;;; Thanks https://stevelosh.com/blog/2021/03/small-common-lisp-cli-programs/
(eval-when (:compile-toplevel :load-toplevel :execute)
(ql:quickload '(:adopt :alexandria :magicl) :silent t))
(defpackage :gather-data
(:use :cl)
(:export :toplevel :*ui*))
(in-package :gather-data)
(defconstant +dtype+ '(complex double-float))
;;;; Modified copy-pasta from quilc
;; https://github.com/quil-lang/quilc/blob/master/src/frontend-utilities.lisp#L161
(defconstant +double-comparison-threshold-loose+ 1d-5)
(defconstant +double-comparison-threshold-strict+ 5d-11)
(defun double~ (x y)
(let ((diff (abs (- x y))))
(< diff +double-comparison-threshold-loose+)))
(defun double= (x y)
(let ((diff (abs (- x y))))
(< diff +double-comparison-threshold-strict+)))
;; https://github.com/quil-lang/quilc/blob/master/src/matrix-operations.lisp#L143
(defun orthonormalize-matrix! (x)
(declare (optimize speed))
(let ((nrows (magicl:nrows x))
(ncols (magicl:ncols x)))
(declare (type alexandria:array-length nrows ncols))
(when (or (zerop nrows) (zerop ncols))
(return-from orthonormalize-matrix! x))
(labels ((column-norm (m col)
(let ((n^2 (dot-columns m col m col)))
(assert (not (double= 0.0d0 n^2)))
(sqrt n^2)))
(normalize-column! (m col)
(let ((norm (column-norm m col)))
(dotimes (row nrows)
(setf (magicl:tref m row col) (/ (magicl:tref m row col) norm)))))
(assign-column! (a b col)
;; a[:,col] = b[:,col]
(dotimes (i nrows)
(setf (magicl:tref a i col) (magicl:tref b i col))))
(dot-columns (a p b q)
(loop :for i :below nrows
:sum (* (conjugate (magicl:tref a i p))
(magicl:tref b i q)))))
(let ((v (magicl:deep-copy-tensor x))
(q (magicl:zeros (list nrows ncols) :type (magicl:element-type x))))
(loop :for j :below ncols :do
(assign-column! q v j)
(normalize-column! q j)
(loop :for k :from (1+ j) :below ncols :do
(let ((s (dot-columns q j v k)))
(loop :for row :below nrows :do
(decf (magicl:tref v row k)
(* s (magicl:tref q row j)))))))
q))))
;; https://github.com/quil-lang/quilc/blob/master/src/compilers/approx.lisp#L158
(defun ensure-positive-determinant (m)
(let ((d (magicl:det m)))
(if (double= -1d0 (realpart d))
(magicl:@ m (magicl:from-diag (list -1 1 1 1) :type +dtype+))
m)))
(defconstant +diagonalizer-max-attempts+ 16)
(defun diagonalizer-number-generator (k)
(abs (sin (/ k 10.0d0))))
(defun orthogonal-diagonalization (uut)
(loop :for attempt :from 0 :below +diagonalizer-max-attempts+ :do
(let* ((coeff (diagonalizer-number-generator attempt))
(matrix (magicl:map (lambda (z)
(+ (* coeff (realpart z))
(* (- 1 coeff) (imagpart z))))
uut))
(evecs (ensure-positive-determinant
(orthonormalize-matrix!
(nth-value 1 (magicl:eig matrix)))))
(evals (magicl:from-diag (magicl:diag (magicl:@ (magicl:transpose evecs) uut evecs))))
(reconstructed (magicl:@ evecs evals (magicl:transpose evecs))))
(when (and (double= 1.0d0 (magicl:det evecs))
(magicl:every #'double= uut reconstructed)
(magicl:every
#'double~
(magicl:eye 4 :type +dtype+)
(magicl:@ (magicl:transpose evecs) evecs)))
(return-from
orthogonal-diagonalization
(list attempt evals evecs)))))
(list +diagonalizer-max-attempts+))
;;;; Functionality
(defconstant +nan+ "NAN")
(defun max-abs-diff (a b)
(loop :for z :across (magicl::storage (magicl:.- a b)) :maximize (abs z)))
(defun collect (u)
(let* ((uut (magicl:@ u (magicl:transpose u)))
(decomp (orthogonal-diagonalization uut)))
(if (< (first decomp) +diagonalizer-max-attempts+)
(destructuring-bind (attempts d o) decomp
(let* ((m (magicl:@ (magicl:transpose o) uut o))
(max-off-diagonal (max-abs-diff
m
(magicl:from-diag (magicl:diag m))))
(max-difference (max-abs-diff
(magicl:@ o d (magicl:transpose o))
uut)))
(list attempts max-off-diagonal max-difference)))
(list +diagonalizer-max-attempts+ +nan+ +nan+))))
(defun report (u name platform stream)
(destructuring-bind (attempts max-off-diagonal max-difference) (collect u)
(format
stream
"~{~A~^,~}~%"
(list
name
attempts
(substitute #\e #\d (format nil "~A" max-off-diagonal))
(substitute #\e #\d (format nil "~A" max-difference))
platform))))
(defconstant
+specifics+
(let* ((i (complex 0d0 1d0))
(-i (complex 0d0 -1d0))
(1/sqrt2 (/ (sqrt 2d0)))
(-i/sqrt2 (* -i 1/sqrt2)))
(list
(magicl:from-list '(1 0 0 0 0 0 0 1 0 0 1 0 0 1 0 0) '(4 4) :type +dtype+)
(magicl:from-diag (list 1 1 1 i) :type +dtype+)
(magicl:kron
(magicl:eye 2 :type +dtype+)
(magicl:from-list (list 1/sqrt2 -i/sqrt2 -i/sqrt2 1/sqrt2) '(2 2) :type +dtype+)))))
(defconstant +iterations+ 1000)
(defun run (platform)
(let ((filename (format nil "~A_lisp.csv" platform)))
(with-open-file (stream filename :direction :output :if-exists :supersede)
(loop :for u :in +specifics+ :do (report u "Specifics" platform stream))
(loop :repeat +iterations+
:for u = (magicl:random-unitary 4)
:then (magicl:random-unitary 4)
:do (report u "SU(4)" platform stream))
(loop :repeat +iterations+
:for u = (magicl:kron (magicl:random-unitary 2) (magicl:random-unitary 2))
:then (magicl:kron (magicl:random-unitary 2) (magicl:random-unitary 2))
:do (report u "SU(2)xSU(2)" platform stream)))))
;;;; CLI
(defparameter *help*
(adopt:make-option 'help
:help "display help and exit"
:long "help"
:short #\h
:reduce (constantly t)))
(defparameter *default-platform* "laptop")
(defparameter *platform*
(adopt:make-option 'platform
:help (format nil "Gather data for PLATFORM (default ~A)" *default-platform*)
:long "platform"
:short #\p
:parameter "PLATFORM"
:initial-value *default-platform*
:reduce #'adopt:last))
(defparameter *ui*
(adopt:make-interface
:name "gather-data"
:usage "[-p PLATFORM]"
:summary "Gather data for the given PLATFORM."
:help "Gather data for the given PLATFORM."
:contents (list *help* *platform*)))
(defun toplevel ()
(handler-case
(multiple-value-bind (arguments options) (adopt:parse-options *ui*)
(when (gethash 'help options)
(adopt:print-help-and-exit *ui*))
(unless (null arguments)
(error "Unrecognized command-line arguments: ~S" arguments))
(run (gethash 'platform options)))
(error (c) (adopt:print-error-and-exit c))))Notebook Diagonalization_Comparison.ipynb (exported to Python) used to generate plots
# |Platform| LAPACK | Installed via |
# |--|--|--|
# |Apple Laptop (Intel) | `liblapack.3.10.1.dylib` | `brew install lapack` |
# | x86 | `liblapack-dev` version `3.7.1-4ubuntu1` | `apt-get install liblapack-dev` |
# | ARM | `liblapack-dev` version `3.9.0-1build1` | `apt-get install liblapack-dev` |
# In[1]:
import matplotlib.pyplot as plt
import pandas as pd
import seaborn as sns
sns.set_theme(
context="notebook",
palette=["#3d3d53", "#00b5ad", "#ef476f"],
style="whitegrid",
rc={"font.sans-serif": "Open Sans"},
)
platforms = ["laptop", "x86", "arm"]
# In[2]:
def plot_attempts(df: pd.DataFrame, title: str):
g = sns.catplot(
df[["name", "attempts", "platform"]],
col="name",
x="attempts",
hue="platform",
kind="count",
sharey=False,
)
g.fig.subplots_adjust(top=0.85)
g.fig.suptitle(title)
def plot_residuals(df: pd.DataFrame, title):
g = sns.catplot(
df.melt(
id_vars=["name", "platform"],
value_vars=[c for c in df.columns if c.startswith("max_")],
),
x="name",
y="value",
col="variable",
hue="platform",
kind="boxen",
)
g.fig.subplots_adjust(top=0.85)
g.fig.suptitle(title)
for ax in g.fig.axes:
ax.set_yscale("log")
# In[3]:
python = pd.concat([pd.read_feather(f"{p}_python.arrow") for p in platforms])
# In[4]:
plot_attempts(python, title="Python Attempts")
# In[5]:
plot_residuals(python, title="Python Residuals")
# In[6]:
plot_residuals(python[python.name != "Specifics"], title="Python Residuals")
# In[7]:
lisp = pd.concat(
[pd.read_csv(f"{p}_lisp.csv", names=python.columns) for p in platforms]
)
# In[8]:
plot_attempts(lisp, title="List Attempts")
# In[9]:
plot_residuals(lisp, title="Lisp Residuals")
# In[10]:
plot_residuals(lisp[lisp.name != "Specifics"], title="Lisp Residuals")
|
I'm sorry this has been such a headache. On the bright side, this is an excellent bug report. |
|
Here's a super duper lil program that reliably triggers the error: from pyquil import get_qc, Program
qc = get_qc("4q-qvm")
# qc = get_qc("Aspen-11", as_qvm=True) # also fails
control = 0
targets = [1, 2, 3]
# targets = [4, 5, 6] # also fails
p = Program(f"""
CNOT {control} {targets[0]}
CNOT {control} {targets[1]}
CNOT {control} {targets[2]}
""")
e = qc.compile(p)Doesn't seem to be related the target architecture (e.g. Aspen or fully connected), or to the target qubits. The only restriction seems to be that the control qubit be the same. |
|
@notmgsk my hero! Running on the ARM box with It also makes a number of complaints (16 to be exact) about |
Spin1Half
added a commit
to Spin1Half/quilc
that referenced
this issue
Sep 30, 2022
It seems that LAPACK becomes unstable on ARM when attempting to find the eigenvectors and eigenvalues of a diagonal matrix. This instability is triggered in find-diagonalizer-in-e-basis. The operation becomes trivial when run on a diagonal matrix, so we check for this and return the identity to avoid making a potentially unstable LAPACK call. Fixes quil-lang#842.
|
@genos @notmgsk Today @stylewarning and I messed around with this. We did not dig down to find the root cause, but we do have a workaround. LAPACK on ARM seems to only get mad over diagonal matrices, which is a trivial case for find-diagonalizer-in-e-basis, so we just check for diagonal matrices and avoid a LAPACK call if this is the case. |
stylewarning
pushed a commit
that referenced
this issue
Sep 30, 2022
It seems that LAPACK becomes unstable on ARM when attempting to find the eigenvectors and eigenvalues of a diagonal matrix. This instability is triggered in find-diagonalizer-in-e-basis. The operation becomes trivial when run on a diagonal matrix, so we check for this and return the identity to avoid making a potentially unstable LAPACK call. Fixes #842.
|
Thanks so much, @Spin1Half! I'm kind of embarrassed that the fix was so simple but exceedingly grateful for it. |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Hello,
quilcwizards! I've caused an interesting crash on ARM64.In order to get
quilcto compile on ARM64, I had to hack oncffia little bit, so perhaps I've screwed something up there? Anyway, when trying to compile "random-ish but architecture-respecting" circuit à la quanvolutional neural networks for running on either of Rigetti's current QPU offerings, I run into "a violent error" (see below). I think it's raised infind-diagonalizer-in-e-basisviacompress-instructions-in-context, but I bow to your knowledge and expertise.Thanks!
The text was updated successfully, but these errors were encountered: