# Documentation for the density mode of MadGraph
Work In Progress

= Density mode
== Authors
V. Durupt
== Short description
Automatic computation of density matrices and quantum information observables. 
This page contains the technical description of the Python library Density_functions.py used for analysis.

== How to compute density matrices
This module of MadGraph5 allows to compute density matrices automatically in a user-friendly way. After having generated your process
{{{
generate YOUR_PROCESS
output
launch
}}}

you will be able to chose the type of run. To compute density matrices one must chose
{{{
reweight = density
}}}

The user is then asked to edit the different cards needed for the run. The information for the density matrices computation is given in the reweight_card. Examples for different benchmark processes are given in https://arxiv.org/abs/2510.17730. At the end of the run, the density matrix for each event is written in the Les Houches Event (LHE) file.

== How to read the density matrices in the LHE file
The LHE file can be read with the parser lhe_parser.py available in MadGraph. An example of code that would parse the density matrix value for each event is:

{{{#!python
import sys
sys.path.append('PATH TO MADGRAPH')
import madgraph.various.Density_functions as dens
import madgraph.various.lhe_parser as lhe_parser

lhe_path = os.path.join('PATH TO MADGRAPH', Data_Path)
for event in lhe_parser.EventFile(lhe_path):
    density = event.density
    square_density = dens.square_matrix(density)

}}}

With this code, for each event of the file, you read the value `density` which contains the independent coefficients of the density matrix. The utility function `square_matrix` allows to write it in the usual matrix form needed for further analysis.\\
For the rest of this page, we will use the same import names as in this block of code.

== Implementation of the quantum information obseravbles
Now that we have the density matrices for each event, we can compute any quantum information observable that we want. The library `Density_functions.py` contains a non-exaustive list of them that will be listed here and that can be called (example of the concurrence here) via.
{{{#!python
concurrence = dens.Get_Concurrence(square_density)
}}}
We will now give a list of the different quantum information observables currently available in the library.
For more physical and more detailed descriptions of the observables, read https://arxiv.org/abs/2510.17730.

=== Purity
Purity represents how mixed a quantum system is.\\ Available for all systems \\It can be computed with the following syntax:
{{{#!python
purity = dens.Get_Purity(square_density)
}}}
Normalised purity (which only differs in the normalisation factor) can be computed with the following syntax:
{{{#!python
purity_normalised = dens.Get_Normalised_Purity(square_density)
}}}

=== Concurrence
Concurrence is an entanglement marker. Concurrence = 0 means the system is separable, and if Concurrence > 0, the bigger it is the more entangled is the system.\\
Available for qubit-qubit (2 x 2) systems only.\\
It can be computed with the following syntax:
{{{#!python
concurrence = dens.Get_Concurrence(square_density)
}}}

For the qutrit-qutrit (3 x 3) system, lower and upper bounds (the square of these bounds are computed more precisely) are known and be computer with the following syntax:
{{{#!python
lower_bound_concurrence_squared = dens.ConcLB2(square_density, pdg_pos)
upper_bound_concurrence_squared = dens.ConcUB2(square_density, pdg_pos)
}}}
Note that `pdg_pos` is a the list of the `PDG` codes of the particles in the density matrix. If you study a system W^+^ W^-^, we would have `pdg_pos = [24, -24]`


=== Bell test
For the case of qubit-qubit systems, the Bell inequality is very simple to compute.\\
Available for qubit-qubit (2 x 2) systems only.\\
It can be computed with the following syntax:
{{{#!python
bell_test, flag_test = dens.Get_Bell_Test(CTC)
}}}
where `CTC` is the product of the transposed spin correlation matrix with the spin correlation matrix and `flag_test` is a boolean being `True` if `\lambda_1 + \lambda_2 > 1`.

=== D coefficients
The D coefficients are conditions for entanglement that are commonly used for the process `p p > t t~`. They are useful to determinate which singlet/triplet state corresponds to a region of high entanglement.\\
Available for qubit-qubit (2 x 2) systems with zero average polarisation.\\
They can be computed with the following syntax:
{{{#!python
D1, Dn, Dr, Dk, boolD = dens.Get_Dcoef(C)
}}}
where `C` is the spin-correlation matrix and `D1, Dn, Dr, Dk` are the D coefficients and `boolD` is `True` if the system shows entanglement.

=== Entanglement of formation
Entanglement of formation is an entanglement marker that has interesting interpretations in term of von Neumann entropy.
Available for qubit-qubit (2 x 2) systems.\\
They can be computed with the following syntax:
{{{#!python
Ef = dens.Get_Entanglement_Formation(square_density)
}}}

=== Magic
Magic is not an entanglement marker, it quantifies another purely quantum feature of a quantum state: the "non-stabiliserness". \\
Available for qubit-qubit (2 x 2) systems.\\
They can be computed with the following syntax:
{{{#!python
Magic = dens.Magic_Mixed(square_density, n)
}}}
where `n` is the number of particles in the density matrix, this has only been tested for `n = 2`.

=== Negativity
Negativity is an entanglement marker based on measuring how negative are the eigenvalues of the partially transposed density matrix. It is based on the Peres-Horodecki criterion which is necessary and sufficient for qubit-qubit (2 x 2) and qubit-qutrit (2 x 3) systems only. Like purity there are two common normalisations, negativity and logarithmic negativity.\\
Computable for every system but interpretation available only for qubit-qubit (2 x 2) and qubit-qutrit (2 x 3) systems.\\
It can be computed with the following syntax:
{{{#!python
Negativity, LogNegativity = dens.Negativity(square_density, pdg_pos)
}}}
where `pdg_pos` is the list of the `PDG` codes in the density matrix.

=== Mana
Mana is an approximation of magic in systems for which magic is not defined.\\
Computable for d x d (d odd) systems like qutrit-qutrit (3 x 3).\\
It can be computed with the following syntax:
{{{#!python
Mana = dens.Get_Mana(square_density, d1, d2)
}}}
where `d1` and `d2` are the dimensions of the Hilbert space on which each of the two particles in the density matrix are defined. Usually we have `d1 = 3 = d2`.

=== Peres-Horodecki criterion
The Peres-Horodecki critetion is the main theorem we have to determine if a system is separable or entangled. The criterion is necessary but not sufficient for all systems and necessary and sufficient for 2 x 2 and 2 x 3 systems. It tells in these cases that: if any of the eigenvalues of the partially transposed density matrix are negative, then the system is entangled.\\
What we compute is the list of the eigenvalues of the partially transposed density matrix.\\
It can be computed with the following syntax:
{{{#!python
flag_entanglement, eigvals = dens.PeresHorodecki_criterion(square_density, pdg_code)
}}}
where `flag_entanglement` is `True` if the system is entangled and `eigvals` are the eigenvalues.

=== Trace and fidelity distance
To compare two different density matrices (from different models, different parameters, different processes, etc.) we can compute "distances". Their interpretation is that the distance is zero if the two matrices represent the same physical process and the distance grows as the matrices get more and more different. We have two computable distances, the trace distance and the fidelity distance.\\
They can be computed for any system with the following syntax:
{{{#!python
trace_distance = dens.trace_distance(square_matrix1, square_matrix2)
fidelity_distance = dens.Fidelity_distance(square_matrix1, square_matrix2)
}}}