-
Notifications
You must be signed in to change notification settings - Fork 11
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Add qexpand functionality and cleanup scalings #59
Conversation
Codecov ReportAttention: Patch coverage is
Additional details and impacted files@@ Coverage Diff @@
## main #59 +/- ##
==========================================
+ Coverage 66.25% 70.10% +3.85%
==========================================
Files 15 17 +2
Lines 646 649 +3
==========================================
+ Hits 428 455 +27
+ Misses 218 194 -24 ☔ View full report in Codecov by Sentry. |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This is shaping really nicely.
I have a couple of questions, and one additional refactoring request. I am not 100% sure about the refactoring request I have (simplifying the constructors and moving the logic to +
, *
, and ⊗
) so do not hesitate to oppose it.
terms = flattenop(+,[c*obj for (c,obj) in d]) | ||
length(d)==1 ? first(terms) : SAdd{S}(d,Set(terms),terms) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This looks a bit suspicious. There is no one-to-one correspondence between the dictionary, set, and list forms anymore. I think what you want is actually countmap_flatten
as seen here
isempty(nonzero_terms) ? f : SAdd{T}(countmap_flatten(nonzero_terms, SScaled{T})) |
Given that there is already flattening done in the +
operation, I think we can just not do any flattening here. Just take the dictionary and do not worry about processing it. Any normal use of SAdd
, usually through the use of +
will already take care of flattening. This is a fairly general statement: fancy processing should be done in calls to +
and *
, not in constructors.
Also, something I did not notice previously, the variable names are misleading here. It should be (obj, coef) in d
not (coef, obj) in d
.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
For this particular case, unlike the other constructors, we probably still want to keep the inner constructor that makes sure to create the set and vector.
@@ -128,7 +128,7 @@ AB | |||
terms | |||
function SMulOperator(terms) | |||
coeff, cleanterms = prefactorscalings(terms) | |||
coeff*new(cleanterms) | |||
coeff*new(flattenop(*,cleanterms)) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I am not sure this is the best order of operations. You might have (c*A) * (c*((c*A)*(c*A)))
-- for this particular example, flattening followed by prefactor extraction would work better. But I am not sure how general my statement is.
A second request: I know this state is kinda inherited already, but now that we are doing more and more fancy processing on construction, I want us to be a bit more disciplined about what goes in the constructor and what goes in the *
operator. May we just remove the custom constructor and move all of this to *
? Any reason that is a bad idea?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I don't see any reason why this would be a bad idea. My original thought process is that using inner constructor methods would be more idiomatic, but perhaps you're right in that it is too fancy for our purposes. It would also make the source code more clear to future developers if we do the processing in the calls to *
, +
, etc.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Thanks! Please proceed with that plan. For more context: one thing that is pestering the back of my mind is how fancy constructors plagued the development of SymPy. At some point, for various reasons related to more advanced simplification rules, it became useful to be able to construct the most raw unsimplified versions of expressions, but all the __init__
methods were already fancy, so a lot of work had to be put in using another set of "pre-init" methods.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Ah, that's interesting. Just out of curiosity, would there be any performance benefits in simplifying the constructors and moving the processing to conversion functions? I can't find any discussion about this in Julia discourse or on GitHub issue pages.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I do not think there would be anything (besides potentially avoiding situation where we happen to write repetitive canonicalization steps). Nothing inherent to the language though (constructors and other functions are treated the same way by the compiler, which is also why we are permitted to return object of type A from the constructor of type B (strictly speaking we are not really constructing type B anymore))
However, this separation of concerns we are discussing above might in some situation make the functions more "type stable" which permits the compiler to make a bunch of optimizations (mostly removing runtime type checks and allocations).
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Thanks, this is good to know.
``` | ||
""" | ||
@withmetadata struct STensor{T<:QObj} <: Symbolic{T} | ||
terms | ||
function STensor{S}(terms) where S | ||
coeff, cleanterms = prefactorscalings(terms) | ||
coeff * new{S}(cleanterms) | ||
coeff * new{S}(flattenop(⊗,cleanterms)) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
same question about order of operations and moving this to ⊗
@@ -202,7 +202,7 @@ julia> commutator(A, A) | |||
op1 | |||
op2 | |||
function SCommutator(o1, o2) | |||
coeff, cleanterms = prefactorscalings([o1 o2], scalar=true) | |||
coeff, cleanterms = prefactorscalings([o1 o2]) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
same question about dropping the custom constructor and moving this to commutator
@@ -233,7 +232,7 @@ julia> anticommutator(A, B) | |||
op1 | |||
op2 | |||
function SAnticommutator(o1, o2) | |||
coeff, cleanterms = prefactorscalings([o1 o2], scalar=true) | |||
coeff, cleanterms = prefactorscalings([o1 o2]) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
same question about dropping the custom constructor and moving this to anticommutator
@@ -92,4 +92,4 @@ symbollabel(x::SZero) = "𝟎" | |||
basis(x::SZero) = nothing | |||
|
|||
Base.show(io::IO, x::SZero) = print(io, symbollabel(x)) | |||
Base.iszero(x::SymQObj) = isa(x, SZero) | |||
Base.iszero(x::Union{SymQObj, Symbolic{Number}, Symbolic{Complex}}) = isa(x, SZero) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I think this is piracy. Why is this needed?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
My bad about the piracy. It was intended for when we process *
for SScaled
objects, although it's kinda of a hackish way to deal with this issue. Basically, I was getting errors before when I would multiply an inner product (which is a subtype of Symbolic{Complex}
) or a real number plus an inner product (subtype of Symbolic{Number}
), as iszero
was not defined for SBraKet
. I'll fix this.
const I = IdentityOp(qubit_basis) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
there seems to be a bit of extra whitespace here
## | ||
|
||
RULES_EXPAND = [ | ||
@rule(commutator(~o1, ~o2) => (~o1)*(~o2) - (~o2)*(~o1)), |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Are the tildas necessary on the right side of the mapping?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Yes. See docs: https://docs.sciml.ai/SymbolicUtils/stable/manual/rewrite/.
@@ -31,10 +30,14 @@ operation(x::SScaled) = * | |||
head(x::SScaled) = :* | |||
children(x::SScaled) = [:*,x.coeff,x.obj] | |||
function Base.:(*)(c, x::Symbolic{T}) where {T<:QObj} | |||
if iszero(c) || iszero(x) | |||
if (isa(c, Number) && iszero(c)) || iszero(x) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This change fixes the type piracy issue with iszero
. Now the *
call checks if the coefficient is a number (and not a symbolic complex number such as an inner product) before it checks if it is zero.
if isa(x, SScaledOperator) | ||
coeff *= x.coeff | ||
push!(terms, x.obj) | ||
elseif isa(x, Union{Number, Symbolic{Number}}) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
In case you were wondering why I added this conditional: when we call flattenop(*, terms)
, we flatten not only multiplicative products, but scalar products as well. So when this is the case, we need to check if there is a scalar value in any of the terms that flattenop
gives. If this wasn't here, then typeof(2*A*B) == SMulOperator
rather than SScaledOperator
. The latter is what we would want, in my opinion.
looks good! Feel free to update the changelog / merge / release (or to postpone some of these steps if it is more convenient for your near term plans) |
This PR addresses several discussions in #55:
prefactorscalings
evolved from its original purpose in a rather messy way. Now it straightforwardly pulls out scalar factors in multiplication and tensor products as well as commutators and anticommutators.flattenop
function that parses repeated addition, multiplication, and tensor products as a single expression. This prevents any ugly nested expressions for these operations.qexpand
function that manually expands expressions containing quantum objects. It is built on top ofSymbolicUtils
utilities for rule-rewriting in a way similar toqsimplify
.utils.jl
andlinalg.jl
and did some rearranging of the source code. For instance, it didn't make sense fordagger
to be inpredefined.jl
where we define things such as Pauli operators and harmonic oscillator operators, sodagger
was moved tolinalg.jl
.express.md
(which was merged in express, qsimplify, and latexify updates #57) to docs.