I've opened a lot of issues and done a lot of prototyping for various ways to improve basis handling here and in dependents. Rather than commenting on each issue separately, I wanted to collect them here and outline my proposal for what a better system looks like and how to get us there. This is partly since there's a lot of subtle interdependence between the issues and fixes for them.

Let's start with the abstract type interface which has been discussed in #26. Based on the experiments in #34 indicating the performance pitfalls of compilation overhead due to combinatorial explosion in type specialization (also see relevant Julia docs), I now think we should not have any type parameters here. If they are necessary in dependent packages, they can reintroduce them easily. So we would have

abstract type Basis end
abstract type StateVector end
abstract type AbstractKet <: StateVector end
abstract type AbstractBra <: StateVector end
abstract type AbstractOperator end
const AbstractQObjType = Union{<:StateVector,<:AbstractOperator} # For convenience
# Note no AbstractSuperOperator, I will address this below

For all the subtypes of Basis defined here, I think we should move to a method interface (see relevant Julia style guide section) which looks like

function basis end # current behavior throwing IncompatibleBases when basis_l(op:T) != basis_r(op:T) for T<:AbstractOperator
function basis_l end # defined for subtypes of AbstractOperator
function basis_r end # defined for subtypes of AbstractOperator
function bases end # defined for subtypes of CompositeBasis and SumBasis and gets the Vector of bases
function spinnumber end # defined for subtypes of SpinBasis
function cutoff end # defined for subtypes of FockBasis
function offset end # defined for subtypes of FockBasis

In addition all subtypes of Basis will be required to implement Base.:(==) and Base.size where the latter must return a tuple with tensor product decomposition of the Hilbert space that the basis defines so that Base.length(b::Basis) = prod(size(b)) gives the total dimension of the Hilbert space (this is currently the .shape field).

Since removing the basis type parameter from AbstractQObjType is meant to go along with moving basis compatibility checks to be entirely dynamic, I think that interface can be improved to look like:

• multiplicable(a::T,b::T) where {T<:AbstractQObjType} checks that they can be multiplied and multiplicable(a,a) checks that an operator is square.
• addible(a::T,b::T) where {T<:AbstractQObjType} checks that they can be added or that the objects can be compared for equality.
• check_* versions of the above with a @compatiblebases macro which disables them.

Finally to address the issue of superoperators and operator bases such as the pauli operator basis (#35, #36). I think I've unfortunately been fairly confused myself about the best way to handle all of this for awhile, but I think I'm hopefully thinking more clearly about it now. I think it's best to view superoperators as simply operators, while the input and output operators they map between are viewed as states corresponding to the vec'd operators.

Concretely this means that the SuperOperator type in QuantumOpticsBase should subtype AbstractOperator and we should introduce a new KetBraBasis <: Basis here. Then we can have that basis_{l,r}(sop::SuperOperator) isa KetBraBasis. More generally all the different types of superoperator representations are related by the different operator bases that they map between. So we can introduce appropriate subtypes of Basis to disambiguate left and right bases of SuperOperator, ChoiState, KrausOperators, PauliTransferMatrix, ChiMatrix. I think this will help a lot with code reuse between the operators*.jl, superoperators.jl and pauli.jl files since there's currently a lot of re-implemented functionality for dense and sparse representations between all those files. I think this also fits much better with Julia's type system. Some operations are the same between operators and superoperators and thus can be shared such as their basis checks. Other are different, such a multiplying two ChoiStates should probably represent the composition of their maps rather than the multiplication of their density matrices even though both are entirely valid things to do with those objects. Method specialization makes this easy to do.

Thoughts @Krastanov @amilsted @apkille?