From 5d2249e1fa9c00dcfdf1831a9efd12703c1390f6 Mon Sep 17 00:00:00 2001 From: Alexander Rosenberg Date: Sat, 7 Dec 2024 21:24:22 -0800 Subject: [PATCH] Continue working on the project in general --- circuit.lisp | 127 +++++++++++++++++++++++++++++++++ math.lisp | 148 ++++++++++++++++++-------------------- package.lisp | 86 +++++++++++++++++++++- parse.lisp | 129 +++++++++++++++++++++++++++++++++ pprint.lisp | 66 +++++++++++++++++ state.lisp | 196 ++++++++++++++++++++------------------------------- 6 files changed, 551 insertions(+), 201 deletions(-) create mode 100644 circuit.lisp create mode 100644 parse.lisp create mode 100644 pprint.lisp diff --git a/circuit.lisp b/circuit.lisp new file mode 100644 index 0000000..c9bcf56 --- /dev/null +++ b/circuit.lisp @@ -0,0 +1,127 @@ +(in-package :cl-quantum/circuit) + +(defun circuit-matrix-operation (state operator &rest args) + (let ((matrix (case operator + (:* identity-2x2) + (:x pauli-x-gate) + (:y pauli-y-gate) + (:z pauli-z-gate) + (:h hadamard-gate) + (:p phase-gate) + (:t pi/8-gate) + (:cnot pauli-x-gate) + (:cz pauli-z-gate)))) + (destructuring-bind (target &optional control) args + (if control + (napply-controlled-operator state matrix target control) + (napply-operator state matrix target))))) + +(defparameter *circuit-operators* + ;; Operator, # args, has output?, function + ;; The output is always the last argument + `((:* 1 nil circuit-matrix-operation) + (:x 1 nil circuit-matrix-operation) + (:y 1 nil circuit-matrix-operation) + (:z 1 nil circuit-matrix-operation) + (:h 1 nil circuit-matrix-operation) + (:p 1 nil circuit-matrix-operation) + (:t 1 nil circuit-matrix-operation) + (:cnot 2 nil circuit-matrix-operation) + (:cz 2 nil circuit-matrix-operation) + (:measure 2 t ,(lambda (state operator &rest args) + (declare (ignorable operator)) + (nmeasure state (car args)))))) + +(defun make-circuit () + "Create a new blank circuit." + '(0)) + +(defun add-to-circuit (circuit operator &rest args) + "Add OPERATOR to CIRCUIT." + (let ((entry (assoc operator *circuit-operators*)) + (largest-arg (apply 'max (remove-if-not 'integerp args)))) + (unless entry + (error "Unknown circuit operator: ~s" operator)) + (destructuring-bind (name arg-count &rest r) entry + (declare (ignorable name r)) + (unless (= arg-count (length args)) + (error "Operator ~s expects ~s args, got ~s" operator + arg-count (length args))) + (when (> (1+ largest-arg) (car circuit)) + (setf (car circuit) (1+ largest-arg))) + (nconc circuit (list (cons operator args)))) + circuit)) + +(defmacro with-circuit (&body body) + "Create a circuit using a simple DSL. BODY can be any valid Lisp forms, in +addition to function calls to functions named in `*circuit-operators*'." + (let ((circuit-var (gensym)) + (size-var (gensym)) + (ao-func (gensym)) + (ao-entry-var (gensym)) + (ao-i-var (gensym))) + `(let ((,circuit-var) + (,size-var 0)) + (flet ((,ao-func (,ao-entry-var) + (dolist (,ao-i-var ,ao-entry-var) + (when (and (integerp ,ao-i-var) + (< ,size-var ,ao-i-var)) + (setq ,size-var ,ao-i-var))) + (push ,ao-entry-var ,circuit-var))) + (macrolet + (,@(mapcar (lambda (oper) + (let ((arg (gensym))) + `(,(car oper) (&rest ,arg) + (assert (= (length ,arg) ,(second oper)) + () + "~s expects ~s arguments, got ~s" + ,(car oper) ,(second oper) (length ,arg)) + `(,',ao-func (list ,',(car oper) ,@,arg))))) + *circuit-operators*)) + ,@body)) + (cons (1+ ,size-var) (nreverse ,circuit-var))))) + +(defun apply-circuit-operator-to-state (state operator args) + "Apply the circuit operator OPERATOR to STATE by calling its function with +ARGS." + (destructuring-bind (&optional name arg-count has-output function) + (assoc operator *circuit-operators*) + (declare (ignorable name arg-count)) + (assert function () + "Invalid circuit operator: ~s" operator) + (let ((output (apply function state operator args))) + (when has-output + (cons (car (last args)) output))))) + +(defun run-circuit (circuit &key bits coefficients probabilities) + "Run the circuit CIRCUIT and return the final state. The initial state can be +specified in one of three ways: + - BITS: the number of qbits + - COEFFICIENTS: the initial coefficients + - PROBABILITES: the initial probabilities" + (assert (= 1 (count-if 'identity (list bits coefficients probabilities))) + () + "Exactly one of BITS, COEFFICIENTS, and PROBABILITIES can be present") + (let ((state (cond + (bits (make-uniform-normal-state bits)) + (coefficients (coerce coefficients 'vector)) + (probabilities (make-normal-state probabilities))))) + (destructuring-bind (circuit-size &rest elements) circuit + (when (> circuit-size (state-bits state)) + (error "Circuit needs at least ~s bits, got ~s" circuit-size + (state-bits state))) + (values + state + (loop for element in elements + for name = (car element) + for args = (cdr element) + for result = (apply-circuit-operator-to-state state name args) + when result + collect result))))) + +(let ((circuit + (with-circuit + (:h 0) + (:cnot 0 1) + (:measure 0 :v1)))) + (run-circuit circuit :bits 2)) diff --git a/math.lisp b/math.lisp index 5c529e2..34ff413 100644 --- a/math.lisp +++ b/math.lisp @@ -1,4 +1,4 @@ -(in-package :cl-quantum) +(in-package :cl-quantum/math) (defmacro domatrix ((var matrix &optional retval) &body body) "Execute BODY for with VAR bound once for each element in MATRIX, then @@ -42,7 +42,7 @@ VALUE, the ROW, and the COLUMN." (setf (aref new-mat row col) (funcall function elem row col))))) ;; Matrix subroutines -(defun mat-minor (mat i j) +(defun minor (mat i j) "Find the minor of MAT for I and J." (destructuring-bind (height width) (array-dimensions mat) @@ -66,7 +66,7 @@ VALUE, the ROW, and the COLUMN." (defun cofactor (mat i j) "Find the cofactor for I and J in MAT." - (* (cofactor-sgn i j) (det (mat-minor mat i j)))) + (* (cofactor-sgn i j) (det (minor mat i j)))) (defun first-column-cofactors (mat) "Find the cofactors for the first column of MAT." @@ -375,8 +375,27 @@ not destructive." (defun round-to-place (num places &key (base 10)) "Round NUM to PLACES places in BASE." - (let ((scale (expt base places))) - (/ (floor (+ (* num scale) 1/2)) scale))) + (if (complexp num) + (let ((real (round-to-place (realpart num) places :base base)) + (imag (round-to-place (imagpart num) places :base base))) + (if (zerop imag) + real + (complex real imag))) + (let ((scale (expt base places))) + (float (/ (floor (+ (* num scale) 1/2)) scale))))) + +(defun round-vector (vec places) + "Round each entry in VEC to PLACES." + (map 'vector (lambda (elt) + (round-to-place elt places)) + vec)) + +(defun round-matrix (mat places) + "Round each entry in MAT to PLACES." + (mapmatrix (lambda (val row col) + (declare (ignorable row col)) + (round-to-place val places)) + mat)) (defun count-digits (num &key (base 10)) "Count the number of digits in NUM. If NUM is zero, return 1. If NUM is @@ -388,80 +407,53 @@ negative, return the number of digits in its absolute value." (defun build-float (int dec) "Create a float with integer part INT and decimal part DEC." - (* (signum int) (+ (abs int) (/ dec (expt 10 (count-digits dec)))))) + (* (if (zerop int) 1 (signum int)) + (+ (abs int) (/ dec (expt 10 (count-digits dec)))))) -(defconstant +parse-real-regexp+ - (ppcre:create-scanner - "^(\\s*([-+]?[0-9]+)(?:/([0-9]+)|\\.?([0-9]*)(?:[eE]([-+]?[0-9]+))?)\\s*)" - :extended-mode t) - "The regexp scanner used in `parse-real'.") +(defun mexp (mat &key (times 100)) + "Calculate exp(MAT) using a Taylor series. The calculation is performed TIMES +times." + (assert (squarep mat) + (mat) + "Matrix must be square: ~s" mat) + (loop for i from 0 to times + for numer = (make-identity-matrix (array-dimension mat 0)) + then (*mm numer mat) + for denom = 1 then (* denom i) + for res = (/ms numer denom) then (+mm res (/ms numer denom)) + finally (return res))) -(defun parse-real (string &key (start 0) end junk-allowed) - "Parse STRING into a real. Parsing starts at START and ends at END. If end is -nil, the end of the string is used. If JUNK-ALLOWED is non-nil, don't signal an -error if an unexpected character is encountered. Two values are returned, the -first being the value parsed and the second being the index at which parsing -stopped. That is, the index of the first un-parsed character." - (values-list - (or - (ppcre:register-groups-bind (whole main denom decim exp) - (+parse-real-regexp+ string :start start :end end :sharedp t) - (unless (or junk-allowed - (= (length whole) (- (or end (length string)) start))) - (error "Malformed number: ~s" (subseq string start end))) - (let ((num - (cond - (denom - (/ (parse-integer main) - (parse-integer denom))) - ((/= (length decim) 0) - (build-float (parse-integer main) - (parse-integer decim))) - (t - (parse-integer main))))) - (list (if exp - (* num (expt 10 (parse-integer exp))) - num) - (length whole)))) - (if junk-allowed - (list 0 0) - (error "Malformed number: ~s" (subseq string start end)))))) +(defun wholep (num) + "Return non-nil if NUM is a whole number." + (or (integerp num) + (and (zerop (imagpart num)) ;; a complex number is not whole + (zerop (second (multiple-value-list (floor (realpart num)))))))) -(defconstant +parse-complex-regexp+ - (ppcre:create-scanner - "^\\s*([-+])?\\s*([-+]?)([0-9/.]+(?:[eE][-+]?[0-9]+)?)?(i)?" - :extended-mode t) - "The regexp scanner used in `parse-complex'.") +(defun mexpt (mat power) + "Calculate MAT to the POWERth power. POWER must be an integer." + (assert (wholep power) + (power) + "Not a whole number: ~s" power) + (let ((acc (make-identity-matrix (array-dimension mat 0)))) + (dotimes (i (floor power) acc) + (setq acc (*mm acc mat))))) -(defun parse-complex (string &key (start 0) end junk-allowed) - "Parse STRING into a complex number. Parsing starts at START and ends at -END. If end is nil, the end of the string is used. If JUNK-ALLOWED is non-nil, -don't signal an error if an unexpected character is encountered. Two values are -returned, the first being the value parsed and the second being the index at -which parsing stopped. That is, the index of the first un-parsed character." - (unless end (setq end (length string))) - (loop for pos = start then (+ pos (length whole)) - for (whole matches) = (multiple-value-list - (ppcre:scan-to-strings +parse-complex-regexp+ - string - :start pos - :end end)) - for times below 2 - while whole - for coef = (cond - ((aref matches 2) - (parse-real (concatenate 'string (aref matches 1) - (aref matches 2)))) - ((aref matches 3) - (if (equal (aref matches 1) "-") -1 1)) - (t 0)) - for sign = (if (equal (aref matches 0) "-") -1 1) - when (aref matches 3) - summing (complex 0 (* sign coef)) into num - else - summing (* sign coef) into num - finally - (if (and (not junk-allowed) - (< pos end)) - (error "Junk in string: ~s" (subseq string start end)) - (return (values num pos))))) +(defun matrix= (m1 m2) + "Return non-nil if each element of M1 and M2 are equal." + (and (= (array-rank m1) (array-rank m2)) + (loop for d1 in (array-dimensions m1) + for d2 in (array-dimensions m2) + unless (= d1 d2) + do (return nil) + finally (return t)) + (domatrix ((row col) m1 t) + (unless (= (aref m1 row col) + (aref m2 row col)) + (return-from matrix=))))) + +(defun vector= (v1 v2) + "Return non-nil if each element of V1 and V2 are equal." + (and (= (length v1) (length v2)) + (dotimes (i (length v1) t) + (unless (= (aref v1 i) (aref v2 i)) + (return-from vector=))))) diff --git a/package.lisp b/package.lisp index ca8f737..a2560ff 100644 --- a/package.lisp +++ b/package.lisp @@ -1,2 +1,84 @@ -(defpackage :cl-quantum - (:use :cl)) +(defpackage :cl-quantum/math + (:use :cl) + (:export #:domatrix + #:mapmatrix + #:minor + #:cofactor + #:det + #:invert + #:transpose + #:norm + #:*mm + #:*mv + #:*vm + #:dot + #:+mm + #:+vv + #:-mm + #:-vv + #:*ms + #:*vs + #:/ms + #:/vs + #:mconj + #:vconj + #:squarep + #:singularp + #:mtrace + #:make-identity-matrix + #:copy-matrix + #:nswap-rows + #:swap-rows + #:nscale-row + #:scale-row + #:nreplace-row-with-sum + #:replace-row-with-sum + #:tensor-mm + #:tensor-vv + #:round-to-place + #:round-vector + #:round-matrix + #:count-digits + #:build-float + #:mexp + #:wholep + #:mexpt + #:matrix= + #:vector=)) + +(defpackage :cl-quantum/state + (:use :cl :cl-quantum/math) + (:export #:normal-state-p + #:normalize-state + #:nnormalize-state + #:state-bits + #:make-normal-state + #:make-uniform-normal-state + #:bit-unset-index + #:bit-set-index + #:bit-probability + #:nmeasure + #:measure + #:make-operator + #:make-controlled-operator + #:apply-operator + #:napply-operator + #:apply-controlled-operator + #:napply-controlled-operator + #:unset-projector + #:set-projector + #:identity-2x2 + #:pauli-x-gate + #:pauli-y-gate + #:pauli-z-gate + #:hadamard-gate + #:phase-gate + #:pi/8-gate + #:cnot-gate + #:cz-gate + #:swap-gate + #:ccnot-gate)) + +(defpackage :cl-quantum/circuit + (:use :cl :cl-quantum/math :cl-quantum/state) + (:export)) diff --git a/parse.lisp b/parse.lisp new file mode 100644 index 0000000..cfba1db --- /dev/null +++ b/parse.lisp @@ -0,0 +1,129 @@ +(defpackage :cl-quantum/parse + (:use :cl :cl-quantum) + (:export #:parse-real + #:parse-complex + #:parse-state + #:parse-bits-state)) + +(in-package :cl-quantum/parse) + +(defconstant +parse-real-regexp+ + (ppcre:create-scanner + "^(\\s*([-+]?[0-9]+)(?:/([0-9]+)|\\.?([0-9]*)(?:[eE]([-+]?[0-9]+))?)\\s*)" + :extended-mode t) + "The regexp scanner used in `parse-real'.") + +(defun parse-real (string &key (start 0) end junk-allowed) + "Parse STRING into a real. Parsing starts at START and ends at END. If end is +nil, the end of the string is used. If JUNK-ALLOWED is non-nil, don't signal an +error if an unexpected character is encountered. Two values are returned, the +first being the value parsed and the second being the index at which parsing +stopped. That is, the index of the first un-parsed character." + (values-list + (or + (ppcre:register-groups-bind (whole main denom decim exp) + (+parse-real-regexp+ string :start start :end end :sharedp t) + (unless (or junk-allowed + (= (length whole) (- (or end (length string)) start))) + (error "Malformed number: ~s" (subseq string start end))) + (let ((num + (cond + (denom + (/ (parse-integer main) + (parse-integer denom))) + ((/= (length decim) 0) + (build-float (parse-integer main) + (parse-integer decim))) + (t + (parse-integer main))))) + (list (if exp + (* num (expt 10 (parse-integer exp))) + num) + (length whole)))) + (if junk-allowed + (list 0 0) + (error "Malformed number: ~s" (subseq string start end)))))) + +(defconstant +parse-complex-regexp+ + (ppcre:create-scanner + "^\\s*([-+])?\\s*([-+]?)([0-9/.]+(?:[eE][-+]?[0-9]+)?)?(i)?" + :extended-mode t) + "The regexp scanner used in `parse-complex'.") + +(defun parse-complex (string &key (start 0) end junk-allowed) + "Parse STRING into a complex number. Parsing starts at START and ends at +END. If end is nil, the end of the string is used. If JUNK-ALLOWED is non-nil, +don't signal an error if an unexpected character is encountered. Two values are +returned, the first being the value parsed and the second being the index at +which parsing stopped. That is, the index of the first un-parsed character." + (unless end (setq end (length string))) + (loop for pos = start then (+ pos (length whole)) + for (whole matches) = (multiple-value-list + (ppcre:scan-to-strings +parse-complex-regexp+ + string + :start pos + :end end)) + for times below 2 + while whole + for coef = (cond + ((aref matches 2) + (parse-real (concatenate 'string (aref matches 1) + (aref matches 2)))) + ((aref matches 3) + (if (equal (aref matches 1) "-") -1 1)) + (t 0)) + for sign = (if (equal (aref matches 0) "-") -1 1) + when (aref matches 3) + summing (complex 0 (* sign coef)) into num + else + summing (* sign coef) into num + finally + (if (and (not junk-allowed) + (< pos end)) + (error "Junk in string: ~s" (subseq string start end)) + (return (values num pos))))) + +(defconstant +parse-state-regexp+ + (ppcre:create-scanner + "^\\s*([-+])?\\s*(\\()?\\s*([-+0-9ei./]*)\\s*(\\))?\\s*\\|([^>]*)>" + :extended-mode t) + "The regexp scanner used in `parse-state'.") + +(defun parse-bits-state (state) + "A `parse-state' parser that parses its state as a binary string." + (parse-integer state :radix 2)) + +(defun parse-state (str &key (parser 'parse-integer)) + "Try to parse STR into a quantum state. PARSER should be a function of one +argument that will take the string inside each ket and return the index of the +state." + (loop for start = 0 then (+ start (length whole)) + for (whole matches) = (multiple-value-list + (ppcre:scan-to-strings +parse-state-regexp+ + str + :sharedp t + :start start)) + while whole + for coef = (if (zerop (length (aref matches 2))) + 1 + (parse-complex (aref matches 2))) + for index = (funcall parser (aref matches 4)) + unless (eq (not (aref matches 1)) + (not (aref matches 3))) + do (error "Mismatches parenthesis: ~s" whole) + when (and (complexp coef) + (not (aref matches 1))) + do (error "Coefficient without matching state: ~s" whole) + collect (if (equal (aref matches 0) "-") + (* -1 coef) + coef) + into coefs + collect index into indecies + maximizing (1+ index) into state-size + finally + (return + (let ((state (make-array state-size :initial-element 0))) + (loop for index in indecies + for coef in coefs + do (incf (aref state index) coef)) + state)))) diff --git a/pprint.lisp b/pprint.lisp new file mode 100644 index 0000000..7f14498 --- /dev/null +++ b/pprint.lisp @@ -0,0 +1,66 @@ +(in-package :cl-quantum) + +(defun pprint-complex (n &key parens (places 5)) + "Pretty-print the complex (or real, rational, etc.) number N. If PARENS is +non-nil, surround the output with parenthesis if it is multiple terms. PLACES is +the number of places to round each part before printing." + (let* ((real (realpart n)) + (imag (imagpart n)) + ;; Also put parenthesis on fractions to make them easier to read + (has-frac (and parens + (or (and (not (integerp real)) + (rationalp real)) + (and (not (integerp imag)) + (rationalp imag)))))) + (when (and (not (rationalp real)) + (realp real)) + (setq real (round-to-place real places))) + (when (and (not (rationalp imag)) + (realp imag)) + (setq imag (round-to-place imag places))) + (cond + ((and (zerop real) + (zerop imag)) + "0") + ((not (or (zerop real) + (zerop imag))) + (format nil "~@[~*(~]~a ~:[-~;+~] ~ai~@[~*)~]" + parens real (>= imag 0) (abs imag) parens)) + (t + (format nil "~@[~*(~]~a~@[~*i~]~@[~*)~]" + has-frac (if (zerop real) imag real) + (zerop real) has-frac))))) + +(defun pprint-format-bits (index size) + "A state formatter that converts the index to binary and pads it with zeros." + (format nil "~v,,,'0<~b~>" (ceiling (log size 2)) index)) + +(defun pprint-format-linear (index size) + "A state formatter that just returns the index +1 as a string." + (declare (ignorable size)) + (format nil "~d" (1+ index))) + +(defun pprint-state (state &key (formatter 'pprint-format-linear) + (places 5)) + "Pretty-print STATE, a quantum state represented as an array. FORMATTER is a +function which takes the index of the quantum state and the total size of the +state. It should convert these to a printable representation. This +representation will be put inside of a ket after each coefficient. PLACES is the +number of places to which to print the coefficients." + (with-output-to-string (out) + (loop with need-sign = nil + for i below (length state) + for coef = (aref state i) + when (and need-sign (not (zerop coef))) + if (>= (realpart coef) 0) + do (format out " + ") + else + do (format out " - ") + and do (setq coef (* -1 coef)) + end + end + unless (zerop coef) + do (format out "~a|~a>" (pprint-complex coef :parens t + :places places) + (funcall formatter i (length state))) + and do (setq need-sign t)))) diff --git a/state.lisp b/state.lisp index de058d4..37d6b2d 100644 --- a/state.lisp +++ b/state.lisp @@ -1,102 +1,4 @@ -(in-package :cl-quantum) - -(defun pprint-complex (n &key parens) - "Pretty-print the complex (or real, rational, etc.) number N. If PARENS is -non-nil, surround the output with parenthesis if it is multiple terms." - (let* ((real (realpart n)) - (imag (imagpart n)) - ;; Also put parenthesis on fractions to make them easier to read - (has-frac (and parens - (or (and (not (integerp real)) - (rationalp real)) - (and (not (integerp imag)) - (rationalp imag)))))) - (cond - ((zerop n) "0") - ((not (or (zerop real) - (zerop imag))) - (format nil "~@[~*(~]~a ~:[-~;+~] ~ai~@[~*)~]" - parens real (>= imag 0) (abs imag) parens)) - (t - (format nil "~@[~*(~]~a~@[~*i~]~@[~*)~]" - has-frac (if (zerop real) imag real) - (zerop real) has-frac))))) - -(defun pprint-format-bits (index size) - "A state formatter that converts the index to binary and pads it with zeros." - (format nil "~v,,,'0<~b~>" (ceiling (log size 2)) index)) - -(defun pprint-format-linear (index size) - "A state formatter that just returns the index +1 as a string." - (declare (ignorable size)) - (format nil "~d" (1+ index))) - -(defun pprint-state (state &key (formatter 'pprint-format-linear)) - "Pretty-print STATE, a quantum state represented as an array. FORMATTER is a -function which takes the index of the quantum state and the total size of the -state. It should convert these to a printable representation. This -representation will be put inside of a ket after each coefficient." - (with-output-to-string (out) - (loop with need-sign = nil - for i below (length state) - for coef = (aref state i) - when (and need-sign (not (zerop coef))) - if (>= (realpart coef) 0) - do (format out " + ") - else - do (format out " - ") - and do (setq coef (* -1 coef)) - end - end - unless (zerop coef) - do (format out "~a|~a>" (pprint-complex coef :parens t) - (funcall formatter i (length state))) - and do (setq need-sign t)))) - -(defconstant +parse-state-regexp+ - (ppcre:create-scanner - "^\\s*([-+])?\\s*(\\()?\\s*([-+0-9ei./]*)\\s*(\\))?\\s*\\|([^>]*)>" - :extended-mode t) - "The regexp scanner used in `parse-state'.") - -(defun parse-bits-state (state) - "A `parse-state' parser that parses its state as a binary string." - (parse-integer state :radix 2)) - -(defun parse-state (str &key (parser 'parse-integer)) - "Try to parse STR into a quantum state. PARSER should be a function of one -argument that will take the string inside each ket and return the index of the -state." - (loop for start = 0 then (+ start (length whole)) - for (whole matches) = (multiple-value-list - (ppcre:scan-to-strings +parse-state-regexp+ - str - :sharedp t - :start start)) - while whole - for coef = (if (zerop (length (aref matches 2))) - 1 - (parse-complex (aref matches 2))) - for index = (funcall parser (aref matches 4)) - unless (eq (not (aref matches 1)) - (not (aref matches 3))) - do (error "Mismatches parenthesis: ~s" whole) - when (and (complexp coef) - (not (aref matches 1))) - do (error "Coefficient without matching state: ~s" whole) - collect (if (equal (aref matches 0) "-") - (* -1 coef) - coef) - into coefs - collect index into indecies - maximizing (1+ index) into state-size - finally - (return - (let ((state (make-array state-size))) - (loop for index in indecies - for coef in coefs - do (incf (aref state index) coef)) - state)))) +(in-package :cl-quantum/state) (defun normal-state-p (state &key (places 5)) "Return non-nil if state is normalized. PLACES is the number of places to @@ -107,6 +9,21 @@ round the norm of STATE before checking." "Return a copy of STATE that is normalized." (/vs state (norm state))) +(defun nnormalize-state (state) + "Normalize STATE by mutating it." + (let ((norm (norm state))) + (dotimes (i (length state) state) + (setf (aref state i) (/ (aref state i) norm))))) + +(defun state-bits (state) + "Return the number of bits in STATE." + (values (ceiling (log (length state) 2)))) + +(defun make-normal-state (probabilities) + "Create a new normalized state with the probability of each state +corresponding to an element in PROBABILITIES." + (map 'vector 'sqrt probabilities)) + (defun make-uniform-normal-state (bits) "Make a uniform normalized quantum state of BITS qbits." (let ((size (ash 1 bits))) @@ -120,7 +37,7 @@ round the norm of STATE before checking." (defun bit-set-index (bit n &key (period (ash 1 bit))) "Return the Nth index in a state in which BIT is 1." - (+ (bit-unset-index bit n :period period))) + (+ period (bit-unset-index bit n :period period))) (defun bit-probability (state bit) "Return the probability that BIT is set in STATE." @@ -131,29 +48,41 @@ round the norm of STATE before checking." for coef = (aref state index) summing (* coef coef))) -(defun nmeasure (state bit) - "Collapse BIT in STATE by measuring it. This will return t or nil depending -on the state the bit collapsed to. Note that this will also modify STATE." - (loop with prob = (round-to-place (bit-probability state bit) 5) +(defun nmeasure (state bit &key (places 5)) + "Collapse BIT in STATE by measuring it. This will return t or nil depending on +the state the bit collapsed to. Note that this will also modify STATE. The +probability will be rounded to PLACES before calculations are carried out." + (loop with prob = (round-to-place (bit-probability state bit) places) + ;; Theocratically there are multiple numbers that could be equivalent + ;; to 100% here, but the chance of those failing is so low anyway that + ;; it doesn't matter. with limit = (* most-positive-fixnum prob) with rnum = (random most-positive-fixnum) - with result = (>= rnum limit) + with result = (< rnum limit) with period = (ash 1 bit) + with 1-prob = (- 1 prob) for i below (/ (length state) 2) for unset-index = (bit-unset-index bit i :period period) for set-index = (+ period unset-index) for unset-coef = (aref state unset-index) for set-coef = (aref state set-index) - for new-coef = (sqrt (+ (* set-coef set-coef) - (* unset-coef unset-coef))) if result do (setf (aref state unset-index) 0 - (aref state set-index) new-coef) + (aref state set-index) + (sqrt (/ (* set-coef set-coef) prob))) else - do (setf (aref state unset-index) new-coef + do (setf (aref state unset-index) + (sqrt (/ (* unset-coef unset-coef) 1-prob)) (aref state set-index) 0) finally (return (values result state)))) +(defun measure (state bit &key (places 5)) + "Like `nmeasure', but don't modify STATE. Note that the new state is returned +as the second value, with the first value being the result of the +measurement. You will probably want to save the second value as the first value +too as the second value alone means pretty much nothing." + (nmeasure (copy-seq state) bit :places places)) + (defun make-operator (bits operator target) "Create an operator matrix that can act on a state with BITS bits and will apply OPERATOR to TARGET." @@ -166,13 +95,13 @@ apply OPERATOR to TARGET." identity-2x2))) finally (return out))) -(defun make-controlled-operator (bits operator target controls) +(defun make-controlled-operator (bits operator target control) "Create an operator matrix that can act on a state with BITS bits and will -apply OPERATOR to TARGET if CONTROLS are all set." +apply OPERATOR to TARGET if CONTROL is set." (labels ((matrix-for (bit target-operator control-operator) (cond ((= bit target) target-operator) - ((member bit controls :test '=) control-operator) + ((= bit control) control-operator) (t identity-2x2))) (tensor-chain (target-operator control-operator) (loop with out = (matrix-for (1- bits) target-operator @@ -184,6 +113,32 @@ apply OPERATOR to TARGET if CONTROLS are all set." (+mm (tensor-chain identity-2x2 unset-projector) (tensor-chain operator set-projector)))) +(defun replace-state (target template) + "Replace each element of TARGET with the corresponding element in TEMPLATE." + (assert (= (length target) (length template)) + (target template) + "Vectors must be of the same length: ~s and ~s" target template) + (dotimes (i (length target) target) + (setf (aref target i) (aref template i)))) + +(defun apply-operator (state operator target) + "Apply OPERATOR to the bit numbered TARGET in STATE." + (*mv (make-operator (state-bits state) operator target) state)) + +(defun napply-operator (state operator target) + "Apply OPERATOR to the bit numbered TARGET in STATE. This modifies state." + (replace-state state (apply-operator state operator target))) + +(defun apply-controlled-operator (state operator target control) + "Apply OPERATOR to the bit numbered TARGET in STATE if CONTROL is set." + (*mv (make-controlled-operator (state-bits state) operator target control) + state)) + +(defun napply-controlled-operator (state operator target control) + "Apply OPERATOR to the bit numbered TARGET in STATE if CONTROL is set. This +modified STATE." + (replace-state state (apply-controlled-operator state operator target control))) + ;;; Gates and Operators: (defconstant unset-projector #2A((1 0) @@ -218,17 +173,16 @@ apply OPERATOR to TARGET if CONTROLS are all set." #2A((1 0) (0 #C(0 1)))) +(defconstant pi/8-gate + (make-array '(2 2) :initial-contents + `((1 0) + (0 ,(exp (complex 0 (/ pi 4))))))) + (defconstant cnot-gate - #2A((1 0 0 0) - (0 1 0 0) - (0 0 0 1) - (0 0 1 0))) + (make-controlled-operator 2 pauli-x-gate 0 1)) (defconstant cz-gate - #2A((1 0 0 0) - (0 1 0 0) - (0 0 1 0) - (0 0 0 -1))) + (make-controlled-operator 2 pauli-z-gate 0 1)) (defconstant swap-gate #2A((1 0 0 0) @@ -237,4 +191,4 @@ apply OPERATOR to TARGET if CONTROLS are all set." (0 0 0 1))) (defconstant ccnot-gate - (nswap-rows (make-identity-matrix 9) 7 8)) + (nswap-rows (make-identity-matrix 8) 6 7))