PIC SETUP:

big box (for 'all' sites/particles) total vol = V length of each side = L V = L³ N_sites = total no.of sites or particles N = total no.of neutrinos at all sites n = total no.density of neutrinos at all sites n = N / V

small box (for 'each' site/particle) total vol = Vᵢ length of each side = Δx Δx³ = Vᵢ Nᵢ = total no.of neutrinos at site i nᵢ = total no.density of neutrinos at site i Nᵢ = nᵢ * Vᵢ

Combining (where index i represent a site and runs from 1:N_sites) ∑ᵢ nᵢ = n # total no.density of neutrinos ∑ᵢ Nᵢ = N # total no.of neutrinos ∑ᵢ Vᵢ = V # total volume of the grid ∑ᵢ Δxᵢ = L # domain size

Create a list of gates that need to be applied to the simulation state.
params: CCNO.Parameters object defined in physics/Parameters.jl
state: CCNO.SimulationSate object defined in physics/SimulationState.jl

Return values are lists of dimensionless unitary operators that include vacuum and self-interactions
gates_1site: list of gates applied to only one site, amenable to simple gate application method
gates_2site_even: list of two-site gates applied to adjacent even pairs of sites. Amenable to simple two-site method.
gates_2site_odd: list of two-site gates applied to adjacent odd pairs of sites. Amenable to simple two-site method.
gates_2site_other: list of two site gates for sites that are not adjacent. Require the general "apply()" function.

Takes the simulation state and advances forward in time according to the parameters passed.

params = CCNO.Parameters object defined in physics/Parameters.jl
state = CCNO.SimulationState object defined in physics/SimulationState.jl

No return value

This function takes in a momentum vector (for all sites) along with their directions 
and returns the geometric factor dependent on the unit momentum vectors for all
i< j sites for a (non) zero p vector.

params: CCNO.Parameters object defined in physics/Parameters.jl
phat: [Nsites,3] array of momentum unit vectors, one unit vector for each site
i,j: indices of the two interacting sites

Returns: scalar geometric factor for the two sites

This function takes in a scalar x  
and returns its heaviside step-function result

This function takes in a momentum vector (for all sites) along with their directions and returns its modulus (for all sites) and unit vectors (for all sites and directions)

Expected (CGS) units of the quantities defined in the files in tests directory that are being used in the gates function.                                                                   
s = site index array (dimensionless and unitless)          
N = array of no.of neutrinos contained on each site (dimensionless and unitless)
B = normalized vector that allows perurbation in different directions (dimensionless constant)
N_sites = Total no.of sites (dimensionless and unitless)
τ = time step (sec)
energy_sign = array of sign of the energy (1 or -1): 1 for neutrinos and -1 for anti-neutrinos
maxdim = max bond dimension in MPS truncation (unitless and dimensionless)
cutoff = truncation threshold for the SVD in MPS representation (unitless and dimensionless)

Return normalized interaction strength used to quantify inhomogeneity.

params: CCNO.Parameters object defined in physics/Parameters.jl
state: CCNO.SimulationState object defined in physics/SimulationState.jl
d: [1,2,3] direction along which the shape function is to be evaluated (x/y/z)
i,j: integers of the two interacting sites