WIP
Some checks failed
MetagraphOptimization_CI / test (push) Failing after 1m29s
MetagraphOptimization_CI / docs (push) Failing after 1m29s

This commit is contained in:
Rubydragon 2024-06-27 13:21:53 +02:00
parent d888713e97
commit 901944bd8b
8 changed files with 710 additions and 337 deletions

View File

@ -14,6 +14,7 @@ JuliaFormatter = "98e50ef6-434e-11e9-1051-2b60c6c9e899"
KernelAbstractions = "63c18a36-062a-441e-b654-da1e3ab1ce7c"
NumaAllocators = "21436f30-1b4a-4f08-87af-e26101bb5379"
QEDbase = "10e22c08-3ccb-4172-bfcf-7d7aa3d04d93"
QEDcore = "35dc0263-cb5f-4c33-a114-1d7f54ab753e"
QEDprocesses = "46de9c38-1bb3-4547-a1ec-da24d767fdad"
Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
Roots = "f2b01f46-fcfa-551c-844a-d8ac1e96c665"

123
examples/congruent_in_ph.jl Normal file
View File

@ -0,0 +1,123 @@
using MetagraphOptimization
using QEDbase
using QEDcore
using Random
using UUIDs
RNG = Random.MersenneTwister(123)
function mock_machine()
return Machine(
[
MetagraphOptimization.NumaNode(
0,
1,
MetagraphOptimization.default_strategy(MetagraphOptimization.NumaNode),
-1.0,
UUIDs.uuid1(),
),
],
[-1.0;;],
)
end
function congruent_input(processDescription::QEDProcessDescription, omega::Number)
# generate an input sample for given e + nk -> e' + k' process, where the nk are equal
massSum = 0
inputMasses = Vector{Float64}()
for (particle, n) in processDescription.inParticles
for _ in 1:n
massSum += mass(particle)
push!(inputMasses, mass(particle))
end
end
outputMasses = Vector{Float64}()
for (particle, n) in processDescription.outParticles
for _ in 1:n
massSum += mass(particle)
push!(outputMasses, mass(particle))
end
end
initialMomenta = [
i == 1 ? SFourMomentum(1, 0, 0, 0) : SFourMomentum(omega, 0, 0, omega) for
i in 1:length(inputMasses)
]
# add some extra random mass to allow for some momentum
ss = sqrt(sum(initialMomenta) * sum(initialMomenta))
result = Vector{QEDProcessInput}()
sizehint!(result, 16)
spin_pol_combinations = Iterators.product(
[SpinUp, SpinDown], [SpinUp, SpinDown], [PolX, PolY], [PolX, PolY]
)
for (in_spin, out_spin, in_pol, out_pol) in spin_pol_combinations
# get the electron first, then the n photons
particles = Vector{QEDParticle}()
for (particle, n) in processDescription.inParticles
if particle <: FermionStateful
mom = initialMomenta[1]
push!(particles, particle(mom, in_spin()))
elseif particle <: PhotonStateful
for i in 1:n
mom = initialMomenta[i + 1]
push!(particles, particle(mom, in_pol()))
end
else
@assert false
end
end
final_momenta = MetagraphOptimization.generate_physical_massive_moms(
RNG, ss, outputMasses
)
index = 1
for (particle, n) in processDescription.outParticles
for _ in 1:n
if particle <: FermionStateful
push!(particles, particle(final_momenta[index], out_spin()))
elseif particle <: PhotonStateful
push!(particles, particle(final_momenta[index], out_pol()))
end
index += 1
end
end
inFerms = MetagraphOptimization._svector_from_type(
processDescription, FermionStateful{Incoming,in_spin}, particles
)
outFerms = MetagraphOptimization._svector_from_type(
processDescription, FermionStateful{Outgoing,out_spin}, particles
)
inAntiferms = MetagraphOptimization._svector_from_type(
processDescription, AntiFermionStateful{Incoming,in_spin}, particles
)
outAntiferms = MetagraphOptimization._svector_from_type(
processDescription, AntiFermionStateful{Outgoing,out_spin}, particles
)
inPhotons = MetagraphOptimization._svector_from_type(
processDescription, PhotonStateful{Incoming,in_pol}, particles
)
outPhotons = MetagraphOptimization._svector_from_type(
processDescription, PhotonStateful{Outgoing,out_pol}, particles
)
processInput = QEDProcessInput(
processDescription,
inFerms,
outFerms,
inAntiferms,
outAntiferms,
inPhotons,
outPhotons,
)
push!(result, processInput)
end
return result
end

View File

@ -6,6 +6,8 @@ A module containing tools to work on DAGs.
module MetagraphOptimization
using QEDbase
using QEDcore
using QEDprocesses
# graph types
export DAG
@ -114,7 +116,6 @@ import Base.delete!
import Base.insert!
import Base.collect
include("devices/interface.jl")
include("task/type.jl")
include("node/type.jl")

View File

@ -5,9 +5,9 @@ using StaticArrays
Return the particle as is and initialize the Value.
"""
function compute(::ComputeTaskQED_P, data::QEDParticleValue{P}) where {P <: QEDParticle}
function compute(::ComputeTaskQED_P, data::QEDParticleValue{P}) where {P<:QEDParticle}
# TODO do we actually need this for anything?
return ParticleValue{P, DiracMatrix}(data.p, one(DiracMatrix))
return ParticleValue{P,DiracMatrix}(data.p, one(DiracMatrix))
end
"""
@ -15,10 +15,12 @@ end
Compute an outer edge. Return the particle value with the same particle and the value multiplied by an outer_edge factor.
"""
function compute(::ComputeTaskQED_U, data::PV) where {P <: QEDParticle, PV <: QEDParticleValue{P}}
function compute(
::ComputeTaskQED_U, data::PV
) where {P<:QEDParticle,PV<:QEDParticleValue{P}}
part::P = data.p
state = base_state(particle(part), direction(part), momentum(part), spin_or_pol(part))
return ParticleValue{P, typeof(state)}(
return ParticleValue{P,typeof(state)}(
data.p,
state, # will return a SLorentzVector{ComplexF64}, BiSpinor or AdjointBiSpinor
)
@ -30,10 +32,10 @@ end
Compute a vertex. Preserve momentum and particle types (e + gamma->p etc.) to create resulting particle, multiply values together and times a vertex factor.
"""
function compute(
::ComputeTaskQED_V,
data1::PV1,
data2::PV2,
) where {P1 <: QEDParticle, P2 <: QEDParticle, PV1 <: QEDParticleValue{P1}, PV2 <: QEDParticleValue{P2}}
::ComputeTaskQED_V, data1::PV1, data2::PV2
) where {
P1<:QEDParticle,P2<:QEDParticle,PV1<:QEDParticleValue{P1},PV2<:QEDParticleValue{P2}
}
p3 = QED_conserve_momentum(data1.p, data2.p)
P3 = interaction_result(P1, P2)
state = QED_vertex()
@ -48,7 +50,7 @@ function compute(
state = state * data2.v
end
dataOut = ParticleValue{P3, typeof(state)}(P3(momentum(p3)), state)
dataOut = ParticleValue{P3,typeof(state)}(P3(momentum(p3)), state)
return dataOut
end
@ -62,10 +64,11 @@ For valid inputs, both input particles should have the same momenta at this poin
12 FLOP.
"""
function compute(
::ComputeTaskQED_S2,
data1::ParticleValue{P1},
data2::ParticleValue{P2},
) where {P1 <: Union{AntiFermionStateful, FermionStateful}, P2 <: Union{AntiFermionStateful, FermionStateful}}
::ComputeTaskQED_S2, data1::ParticleValue{P1}, data2::ParticleValue{P2}
) where {
P1<:Union{AntiFermionStateful,FermionStateful},
P2<:Union{AntiFermionStateful,FermionStateful},
}
#@assert isapprox(data1.p.momentum, data2.p.momentum, rtol = sqrt(eps()), atol = sqrt(eps())) "$(data1.p.momentum) vs. $(data2.p.momentum)"
inner = QED_inner_edge(propagation_result(P1)(momentum(data1.p)))
@ -79,10 +82,8 @@ function compute(
end
function compute(
::ComputeTaskQED_S2,
data1::ParticleValue{P1},
data2::ParticleValue{P2},
) where {P1 <: PhotonStateful, P2 <: PhotonStateful}
::ComputeTaskQED_S2, data1::ParticleValue{P1}, data2::ParticleValue{P2}
) where {P1<:PhotonStateful,P2<:PhotonStateful}
# TODO: assert that data1 and data2 are opposites
inner = QED_inner_edge(data1.p)
# inner edge is just a scalar, data1 and data2 are photon states that are just Complex numbers here
@ -94,7 +95,7 @@ end
Compute inner edge (1 input particle, 1 output particle).
"""
function compute(::ComputeTaskQED_S1, data::QEDParticleValue{P}) where {P <: QEDParticle}
function compute(::ComputeTaskQED_S1, data::QEDParticleValue{P}) where {P<:QEDParticle}
newP = propagation_result(P)
new_p = newP(momentum(data.p))
# inner edge is just a scalar, can multiply from either side

View File

@ -1,14 +1,28 @@
ComputeTaskQED_Sum() = ComputeTaskQED_Sum(0)
function _svector_from_type(processDescription::QEDProcessDescription, type, particles)
if haskey(in_particles(processDescription), type)
return SVector{in_particles(processDescription)[type], type}(filter(x -> typeof(x) <: type, particles))
function _svector_from_type(
processDescription::QEDProcessDescription, type::Type{T}, particles
) where {DIR<:ParticleDirection,T<:QEDParticle{DIR}}
if DIR <: Incoming
l = 0
for (k, v) in in_particles(processDescription)
if T <: k
l = v
break
end
end
return SVector{l,T}(filter(x -> typeof(x) <: T, particles))
elseif DIR <: Outgoing
l = 0
for (k, v) in out_particles(processDescription)
if T <: k
l = v
break
end
end
return SVector{l,T}(filter(x -> typeof(x) <: T, particles))
end
if haskey(out_particles(processDescription), type)
return SVector{out_particles(processDescription)[type], type}(filter(x -> typeof(x) <: type, particles))
end
return SVector{0, type}()
end
"""
@ -38,7 +52,6 @@ function gen_process_input(processDescription::QEDProcessDescription)
# add some extra random mass to allow for some momentum
massSum += rand(rng[threadid()]) * (length(inputMasses) + length(outputMasses))
particles = Vector{QEDParticle}()
initialMomenta = generate_initial_moms(massSum, inputMasses)
index = 1
@ -59,15 +72,34 @@ function gen_process_input(processDescription::QEDProcessDescription)
end
end
inFerms = _svector_from_type(processDescription, FermionStateful{Incoming, SpinUp}, particles)
outFerms = _svector_from_type(processDescription, FermionStateful{Outgoing, SpinUp}, particles)
inAntiferms = _svector_from_type(processDescription, AntiFermionStateful{Incoming, SpinUp}, particles)
outAntiferms = _svector_from_type(processDescription, AntiFermionStateful{Outgoing, SpinUp}, particles)
inPhotons = _svector_from_type(processDescription, PhotonStateful{Incoming, PolX}, particles)
outPhotons = _svector_from_type(processDescription, PhotonStateful{Outgoing, PolX}, particles)
inFerms = _svector_from_type(
processDescription, FermionStateful{Incoming,SpinUp}, particles
)
outFerms = _svector_from_type(
processDescription, FermionStateful{Outgoing,SpinUp}, particles
)
inAntiferms = _svector_from_type(
processDescription, AntiFermionStateful{Incoming,SpinUp}, particles
)
outAntiferms = _svector_from_type(
processDescription, AntiFermionStateful{Outgoing,SpinUp}, particles
)
inPhotons = _svector_from_type(
processDescription, PhotonStateful{Incoming,PolX}, particles
)
outPhotons = _svector_from_type(
processDescription, PhotonStateful{Outgoing,PolX}, particles
)
processInput =
QEDProcessInput(processDescription, inFerms, outFerms, inAntiferms, outAntiferms, inPhotons, outPhotons)
processInput = QEDProcessInput(
processDescription,
inFerms,
outFerms,
inAntiferms,
outAntiferms,
inPhotons,
outPhotons,
)
return processInput
end
@ -88,9 +120,13 @@ function gen_graph(process_description::QEDProcessDescription)
# TODO: Not all diagram outputs should always be summed at the end, if they differ by fermion exchange they need to be diffed
# Should not matter for n-Photon Compton processes though
sum_node = insert_node!(graph, make_node(ComputeTaskQED_Sum(0)), track = false, invalidate_cache = false)
global_data_out = insert_node!(graph, make_node(DataTask(COMPLEX_SIZE)), track = false, invalidate_cache = false)
insert_edge!(graph, sum_node, global_data_out, track = false, invalidate_cache = false)
sum_node = insert_node!(
graph, make_node(ComputeTaskQED_Sum(0)); track=false, invalidate_cache=false
)
global_data_out = insert_node!(
graph, make_node(DataTask(COMPLEX_SIZE)); track=false, invalidate_cache=false
)
insert_edge!(graph, sum_node, global_data_out; track=false, invalidate_cache=false)
# remember the data out nodes for connection
dataOutNodes = Dict()
@ -99,16 +135,22 @@ function gen_graph(process_description::QEDProcessDescription)
# generate data in and U tasks
data_in = insert_node!(
graph,
make_node(DataTask(PARTICLE_VALUE_SIZE), String(particle)),
track = false,
invalidate_cache = false,
make_node(DataTask(PARTICLE_VALUE_SIZE), String(particle));
track=false,
invalidate_cache=false,
) # read particle data node
compute_u = insert_node!(graph, make_node(ComputeTaskQED_U()), track = false, invalidate_cache = false) # compute U node
data_out =
insert_node!(graph, make_node(DataTask(PARTICLE_VALUE_SIZE)), track = false, invalidate_cache = false) # transfer data out from u (one ABCParticleValue object)
compute_u = insert_node!(
graph, make_node(ComputeTaskQED_U()); track=false, invalidate_cache=false
) # compute U node
data_out = insert_node!(
graph,
make_node(DataTask(PARTICLE_VALUE_SIZE));
track=false,
invalidate_cache=false,
) # transfer data out from u (one ABCParticleValue object)
insert_edge!(graph, data_in, compute_u, track = false, invalidate_cache = false)
insert_edge!(graph, compute_u, data_out, track = false, invalidate_cache = false)
insert_edge!(graph, data_in, compute_u; track=false, invalidate_cache=false)
insert_edge!(graph, compute_u, data_out; track=false, invalidate_cache=false)
# remember the data_out node for future edges
dataOutNodes[String(particle)] = data_out
@ -124,19 +166,30 @@ function gen_graph(process_description::QEDProcessDescription)
data_in1 = dataOutNodes[String(vertex.in1)]
data_in2 = dataOutNodes[String(vertex.in2)]
compute_V = insert_node!(graph, make_node(ComputeTaskQED_V()), track = false, invalidate_cache = false) # compute vertex
compute_V = insert_node!(
graph,
make_node(ComputeTaskQED_V());
track=false,
invalidate_cache=false,
) # compute vertex
insert_edge!(graph, data_in1, compute_V, track = false, invalidate_cache = false)
insert_edge!(graph, data_in2, compute_V, track = false, invalidate_cache = false)
insert_edge!(
graph, data_in1, compute_V; track=false, invalidate_cache=false
)
insert_edge!(
graph, data_in2, compute_V; track=false, invalidate_cache=false
)
data_V_out = insert_node!(
graph,
make_node(DataTask(PARTICLE_VALUE_SIZE)),
track = false,
invalidate_cache = false,
make_node(DataTask(PARTICLE_VALUE_SIZE));
track=false,
invalidate_cache=false,
)
insert_edge!(graph, compute_V, data_V_out, track = false, invalidate_cache = false)
insert_edge!(
graph, compute_V, data_V_out; track=false, invalidate_cache=false
)
if (vertex.out == tie.in1 || vertex.out == tie.in2)
# out particle is part of the tie -> there will be an S2 task with it later, don't make S1 task
@ -145,19 +198,27 @@ function gen_graph(process_description::QEDProcessDescription)
end
# otherwise, add S1 task
compute_S1 =
insert_node!(graph, make_node(ComputeTaskQED_S1()), track = false, invalidate_cache = false) # compute propagator
compute_S1 = insert_node!(
graph,
make_node(ComputeTaskQED_S1());
track=false,
invalidate_cache=false,
) # compute propagator
insert_edge!(graph, data_V_out, compute_S1, track = false, invalidate_cache = false)
insert_edge!(
graph, data_V_out, compute_S1; track=false, invalidate_cache=false
)
data_S1_out = insert_node!(
graph,
make_node(DataTask(PARTICLE_VALUE_SIZE)),
track = false,
invalidate_cache = false,
make_node(DataTask(PARTICLE_VALUE_SIZE));
track=false,
invalidate_cache=false,
)
insert_edge!(graph, compute_S1, data_S1_out, track = false, invalidate_cache = false)
insert_edge!(
graph, compute_S1, data_S1_out; track=false, invalidate_cache=false
)
# overrides potentially different nodes from previous diagrams, which is intentional
dataOutNodes[String(vertex.out)] = data_S1_out
@ -168,16 +229,23 @@ function gen_graph(process_description::QEDProcessDescription)
data_in1 = dataOutNodes[String(tie.in1)]
data_in2 = dataOutNodes[String(tie.in2)]
compute_S2 = insert_node!(graph, make_node(ComputeTaskQED_S2()), track = false, invalidate_cache = false)
compute_S2 = insert_node!(
graph, make_node(ComputeTaskQED_S2()); track=false, invalidate_cache=false
)
data_S2 = insert_node!(graph, make_node(DataTask(PARTICLE_VALUE_SIZE)), track = false, invalidate_cache = false)
data_S2 = insert_node!(
graph,
make_node(DataTask(PARTICLE_VALUE_SIZE));
track=false,
invalidate_cache=false,
)
insert_edge!(graph, data_in1, compute_S2, track = false, invalidate_cache = false)
insert_edge!(graph, data_in2, compute_S2, track = false, invalidate_cache = false)
insert_edge!(graph, data_in1, compute_S2; track=false, invalidate_cache=false)
insert_edge!(graph, data_in2, compute_S2; track=false, invalidate_cache=false)
insert_edge!(graph, compute_S2, data_S2, track = false, invalidate_cache = false)
insert_edge!(graph, compute_S2, data_S2; track=false, invalidate_cache=false)
insert_edge!(graph, data_S2, sum_node, track = false, invalidate_cache = false)
insert_edge!(graph, data_S2, sum_node; track=false, invalidate_cache=false)
add_child!(task(sum_node))
end

View File

@ -45,9 +45,9 @@ The [`FeynmanTie`](@ref) represents the final inner edge of the diagram.
"""
struct FeynmanDiagram
vertices::Vector{Set{FeynmanVertex}}
tie::Ref{Union{FeynmanTie, Missing}}
tie::Ref{Union{FeynmanTie,Missing}}
particles::Vector{FeynmanParticle}
type_ids::Dict{Type, Int64} # lut for number of used ids for a particle type
type_ids::Dict{Type,Int64} # lut for number of used ids for a particle type
end
"""
@ -69,7 +69,7 @@ function FeynmanDiagram(pd::QEDProcessDescription)
push!(parts, FeynmanParticle(type, i))
end
end
ids = Dict{Type, Int64}()
ids = Dict{Type,Int64}()
for t in types(QEDModel())
if (isincoming(t))
ids[t] = get(pd.inParticles, t, 0)
@ -83,13 +83,17 @@ end
function particle_after_tie(p::FeynmanParticle, t::FeynmanTie)
if p == t.in1 || p == t.in2
return FeynmanParticle(FermionStateful{Incoming, SpinUp}, -1) # placeholder particle and id for tied particles
return FeynmanParticle(FermionStateful{Incoming,SpinUp}, -1) # placeholder particle and id for tied particles
end
return p
end
function vertex_after_tie(v::FeynmanVertex, t::FeynmanTie)
return FeynmanVertex(particle_after_tie(v.in1, t), particle_after_tie(v.in2, t), particle_after_tie(v.out, t))
return FeynmanVertex(
particle_after_tie(v.in1, t),
particle_after_tie(v.in2, t),
particle_after_tie(v.out, t),
)
end
function vertex_after_tie(v::FeynmanVertex, t::Missing)
@ -104,7 +108,9 @@ function vertex_set_after_tie(vs::Set{FeynmanVertex}, t::Missing)
return vs
end
function vertex_set_after_tie(vs::Set{FeynmanVertex}, t1::Union{FeynmanTie, Missing}, t2::Union{FeynmanTie, Missing})
function vertex_set_after_tie(
vs::Set{FeynmanVertex}, t1::Union{FeynmanTie,Missing}, t2::Union{FeynmanTie,Missing}
)
return Set{FeynmanVertex}(vertex_after_tie(vertex_after_tie(v, t1), t2) for v in vs)
end
@ -138,7 +144,8 @@ function ==(t1::FeynmanTie, t2::FeynmanTie)
end
function ==(d1::FeynmanDiagram, d2::FeynmanDiagram)
if (!ismissing(d1.tie[]) && ismissing(d2.tie[])) || (ismissing(d1.tie[]) && !ismissing(d2.tie[]))
if (!ismissing(d1.tie[]) && ismissing(d2.tie[])) ||
(ismissing(d1.tie[]) && !ismissing(d2.tie[]))
return false
end
if d1.particles != d2.particles
@ -150,7 +157,8 @@ function ==(d1::FeynmanDiagram, d2::FeynmanDiagram)
# TODO can i prove that this works?
for (v1, v2) in zip(d1.vertices, d2.vertices)
if vertex_set_after_tie(v1, d1.tie[], d2.tie[]) != vertex_set_after_tie(v2, d1.tie[], d2.tie[])
if vertex_set_after_tie(v1, d1.tie[], d2.tie[]) !=
vertex_set_after_tie(v2, d1.tie[], d2.tie[])
return false
end
end
@ -162,8 +170,11 @@ function ==(d1::FeynmanDiagram, d2::FeynmanDiagram)
)=#
end
copy(fd::FeynmanDiagram) =
FeynmanDiagram(deepcopy(fd.vertices), copy(fd.tie[]), deepcopy(fd.particles), copy(fd.type_ids))
function copy(fd::FeynmanDiagram)
return FeynmanDiagram(
deepcopy(fd.vertices), copy(fd.tie[]), deepcopy(fd.particles), copy(fd.type_ids)
)
end
"""
id_for_type(d::FeynmanDiagram, t::Type{<:QEDParticle})
@ -229,7 +240,7 @@ end
Return a vector of the particles after applying the vertices and tie of the diagram up to the given level. If no level is given, apply all. The tie comes last and is its own "level".
"""
function get_particles(fd::FeynmanDiagram, level::Int = -1)
function get_particles(fd::FeynmanDiagram, level::Int=-1)
if level == -1
level = length(fd.vertices) + 1
end
@ -365,8 +376,14 @@ function possible_vertices(fd::FeynmanDiagram)
p1 = particles[i]
p2 = particles[j]
if (caninteract(p1.particle, p2.particle))
interaction_res = propagation_result(interaction_result(p1.particle, p2.particle))
v = FeynmanVertex(p1, p2, FeynmanParticle(interaction_res, id_for_type(fd, interaction_res)))
interaction_res = propagation_result(
interaction_result(p1.particle, p2.particle)
)
v = FeynmanVertex(
p1,
p2,
FeynmanParticle(interaction_res, id_for_type(fd, interaction_res)),
)
#@assert !(v.out in particles) "$v is in $fd"
if !can_apply_vertex(fully_generated_particles, v)
continue
@ -445,12 +462,12 @@ end
Returns true iff the given feynman diagram is an (empty) diagram of a compton process like ke->k^ne
"""
function is_compton(fd::FeynmanDiagram)
return fd.type_ids[FermionStateful{Incoming, SpinUp}] == 1 &&
fd.type_ids[FermionStateful{Outgoing, SpinUp}] == 1 &&
fd.type_ids[AntiFermionStateful{Incoming, SpinUp}] == 0 &&
fd.type_ids[AntiFermionStateful{Outgoing, SpinUp}] == 0 &&
fd.type_ids[PhotonStateful{Incoming, PolX}] >= 1 &&
fd.type_ids[PhotonStateful{Outgoing, PolX}] >= 1
return fd.type_ids[FermionStateful{Incoming}] == 1 &&
fd.type_ids[FermionStateful{Outgoing}] == 1 &&
fd.type_ids[AntiFermionStateful{Incoming}] == 0 &&
fd.type_ids[AntiFermionStateful{Outgoing}] == 0 &&
fd.type_ids[PhotonStateful{Incoming}] >= 1 &&
fd.type_ids[PhotonStateful{Outgoing}] >= 1
end
"""
@ -460,19 +477,19 @@ Helper function for [`gen_compton_diagrams`](@Ref). Generates a single diagram f
"""
function gen_compton_diagram_from_order(order::Vector{Int}, inFerm, outFerm, n::Int, m::Int)
photons = vcat(
[FeynmanParticle(PhotonStateful{Incoming, PolX}, i) for i in 1:n],
[FeynmanParticle(PhotonStateful{Outgoing, PolX}, i) for i in 1:m],
[FeynmanParticle(PhotonStateful{Incoming,PolX}, i) for i in 1:n],
[FeynmanParticle(PhotonStateful{Outgoing,PolX}, i) for i in 1:m],
)
new_diagram = FeynmanDiagram(
[],
missing,
[inFerm, outFerm, photons...],
Dict{Type, Int64}(
FermionStateful{Incoming, SpinUp} => 1,
FermionStateful{Outgoing, SpinUp} => 1,
PhotonStateful{Incoming, PolX} => n,
PhotonStateful{Outgoing, PolX} => m,
Dict{Type,Int64}(
FermionStateful{Incoming,SpinUp} => 1,
FermionStateful{Outgoing,SpinUp} => 1,
PhotonStateful{Incoming,PolX} => n,
PhotonStateful{Outgoing,PolX} => m,
),
)
@ -484,9 +501,9 @@ function gen_compton_diagram_from_order(order::Vector{Int}, inFerm, outFerm, n::
while left_index <= right_index
# left side
v_left = FeynmanVertex(
FeynmanParticle(FermionStateful{Incoming, SpinUp}, iterations),
FeynmanParticle(FermionStateful{Incoming,SpinUp}, iterations),
photons[order[left_index]],
FeynmanParticle(FermionStateful{Incoming, SpinUp}, iterations + 1),
FeynmanParticle(FermionStateful{Incoming,SpinUp}, iterations + 1),
)
left_index += 1
add_vertex!(new_diagram, v_left)
@ -497,9 +514,9 @@ function gen_compton_diagram_from_order(order::Vector{Int}, inFerm, outFerm, n::
# right side
v_right = FeynmanVertex(
FeynmanParticle(FermionStateful{Outgoing, SpinUp}, iterations),
FeynmanParticle(FermionStateful{Outgoing,SpinUp}, iterations),
photons[order[right_index]],
FeynmanParticle(FermionStateful{Outgoing, SpinUp}, iterations + 1),
FeynmanParticle(FermionStateful{Outgoing,SpinUp}, iterations + 1),
)
right_index -= 1
add_vertex!(new_diagram, v_right)
@ -512,27 +529,28 @@ function gen_compton_diagram_from_order(order::Vector{Int}, inFerm, outFerm, n::
return new_diagram
end
"""
gen_compton_diagram_from_order_one_side(order::Vector{Int}, inFerm, outFerm, n::Int, m::Int)
Helper function for [`gen_compton_diagrams`](@Ref). Generates a single diagram for the given order and n input and m output photons.
"""
function gen_compton_diagram_from_order_one_side(order::Vector{Int}, inFerm, outFerm, n::Int, m::Int)
function gen_compton_diagram_from_order_one_side(
order::Vector{Int}, inFerm, outFerm, n::Int, m::Int
)
photons = vcat(
[FeynmanParticle(PhotonStateful{Incoming, PolX}, i) for i in 1:n],
[FeynmanParticle(PhotonStateful{Outgoing, PolX}, i) for i in 1:m],
[FeynmanParticle(PhotonStateful{Incoming,PolX}, i) for i in 1:n],
[FeynmanParticle(PhotonStateful{Outgoing,PolX}, i) for i in 1:m],
)
new_diagram = FeynmanDiagram(
[],
missing,
[inFerm, outFerm, photons...],
Dict{Type, Int64}(
FermionStateful{Incoming, SpinUp} => 1,
FermionStateful{Outgoing, SpinUp} => 1,
PhotonStateful{Incoming, PolX} => n,
PhotonStateful{Outgoing, PolX} => m,
Dict{Type,Int64}(
FermionStateful{Incoming,SpinUp} => 1,
FermionStateful{Outgoing,SpinUp} => 1,
PhotonStateful{Incoming,PolX} => n,
PhotonStateful{Outgoing,PolX} => m,
),
)
@ -544,9 +562,9 @@ function gen_compton_diagram_from_order_one_side(order::Vector{Int}, inFerm, out
while left_index <= right_index
# left side
v_left = FeynmanVertex(
FeynmanParticle(FermionStateful{Incoming, SpinUp}, iterations),
FeynmanParticle(FermionStateful{Incoming,SpinUp}, iterations),
photons[order[left_index]],
FeynmanParticle(FermionStateful{Incoming, SpinUp}, iterations + 1),
FeynmanParticle(FermionStateful{Incoming,SpinUp}, iterations + 1),
)
left_index += 1
add_vertex!(new_diagram, v_left)
@ -559,9 +577,9 @@ function gen_compton_diagram_from_order_one_side(order::Vector{Int}, inFerm, out
if (iterations == 1)
# right side
v_right = FeynmanVertex(
FeynmanParticle(FermionStateful{Outgoing, SpinUp}, iterations),
FeynmanParticle(FermionStateful{Outgoing,SpinUp}, iterations),
photons[order[right_index]],
FeynmanParticle(FermionStateful{Outgoing, SpinUp}, iterations + 1),
FeynmanParticle(FermionStateful{Outgoing,SpinUp}, iterations + 1),
)
right_index -= 1
add_vertex!(new_diagram, v_right)
@ -576,41 +594,45 @@ function gen_compton_diagram_from_order_one_side(order::Vector{Int}, inFerm, out
return new_diagram
end
"""
gen_compton_diagrams(n::Int, m::Int)
Special case diagram generation for Compton processes, i.e., processes of the form k^ne->k^me
"""
function gen_compton_diagrams(n::Int, m::Int)
inFerm = FeynmanParticle(FermionStateful{Incoming, SpinUp}, 1)
outFerm = FeynmanParticle(FermionStateful{Outgoing, SpinUp}, 1)
inFerm = FeynmanParticle(FermionStateful{Incoming,SpinUp}, 1)
outFerm = FeynmanParticle(FermionStateful{Outgoing,SpinUp}, 1)
perms = [permutations([i for i in 1:(n + m)])...]
diagrams = [Vector{FeynmanDiagram}() for i in 1:nthreads()]
@threads for order in perms
push!(diagrams[threadid()], gen_compton_diagram_from_order(order, inFerm, outFerm, n, m))
push!(
diagrams[threadid()],
gen_compton_diagram_from_order(order, inFerm, outFerm, n, m),
)
end
return vcat(diagrams...)
end
"""
gen_compton_diagrams_one_side(n::Int, m::Int)
Special case diagram generation for Compton processes, i.e., processes of the form k^ne->k^me, but generating from one end, yielding larger diagrams
"""
function gen_compton_diagrams_one_side(n::Int, m::Int)
inFerm = FeynmanParticle(FermionStateful{Incoming, SpinUp}, 1)
outFerm = FeynmanParticle(FermionStateful{Outgoing, SpinUp}, 1)
inFerm = FeynmanParticle(FermionStateful{Incoming,SpinUp}, 1)
outFerm = FeynmanParticle(FermionStateful{Outgoing,SpinUp}, 1)
perms = [permutations([i for i in 1:(n + m)])...]
diagrams = [Vector{FeynmanDiagram}() for i in 1:nthreads()]
@threads for order in perms
push!(diagrams[threadid()], gen_compton_diagram_from_order_one_side(order, inFerm, outFerm, n, m))
push!(
diagrams[threadid()],
gen_compton_diagram_from_order_one_side(order, inFerm, outFerm, n, m),
)
end
return vcat(diagrams...)
@ -624,8 +646,7 @@ From a given feynman diagram in its initial state, e.g. when created through the
function gen_diagrams(fd::FeynmanDiagram)
if is_compton(fd)
return gen_compton_diagrams_one_side(
fd.type_ids[PhotonStateful{Incoming, PolX}],
fd.type_ids[PhotonStateful{Outgoing, PolX}],
fd.type_ids[PhotonStateful{Incoming}], fd.type_ids[PhotonStateful{Outgoing}]
)
end

View File

@ -21,7 +21,7 @@ Its template parameter specifies the particle's direction.
The concrete types contain singletons of the types that they are, like `Photon` and `Electron` from QEDbase, and their state descriptions.
"""
abstract type QEDParticle{Direction <: ParticleDirection} <: AbstractParticle end
abstract type QEDParticle{Direction<:ParticleDirection} <: AbstractParticle end
"""
QEDProcessDescription <: AbstractProcessDescription
@ -31,16 +31,16 @@ A description of a process in the QED-Model. Contains the input and output parti
See also: [`in_particles`](@ref), [`out_particles`](@ref), [`parse_process`](@ref)
"""
struct QEDProcessDescription <: AbstractProcessDescription
inParticles::Dict{Type{<:QEDParticle{Incoming}}, Int}
outParticles::Dict{Type{<:QEDParticle{Outgoing}}, Int}
inParticles::Dict{Type{<:QEDParticle{Incoming}},Int}
outParticles::Dict{Type{<:QEDParticle{Outgoing}},Int}
end
QEDParticleValue{ParticleType <: QEDParticle} = Union{
ParticleValue{ParticleType, BiSpinor},
ParticleValue{ParticleType, AdjointBiSpinor},
ParticleValue{ParticleType, DiracMatrix},
ParticleValue{ParticleType, SLorentzVector{Float64}},
ParticleValue{ParticleType, ComplexF64},
QEDParticleValue{ParticleType<:QEDParticle} = Union{
ParticleValue{ParticleType,BiSpinor},
ParticleValue{ParticleType,AdjointBiSpinor},
ParticleValue{ParticleType,DiracMatrix},
ParticleValue{ParticleType,SLorentzVector{Float64}},
ParticleValue{ParticleType,ComplexF64},
}
"""
@ -48,84 +48,158 @@ QEDParticleValue{ParticleType <: QEDParticle} = Union{
A photon of the [`QEDModel`](@ref) with its state.
"""
struct PhotonStateful{Direction <: ParticleDirection, Pol <: AbstractDefinitePolarization} <: QEDParticle{Direction}
struct PhotonStateful{Direction<:ParticleDirection,Pol<:AbstractDefinitePolarization} <:
QEDParticle{Direction}
momentum::SFourMomentum
function PhotonStateful{Direction,Pol}(
mom::SFourMomentum
) where {Direction<:ParticleDirection,Pol<:AbstractDefinitePolarization}
return new{Direction,Pol}(mom)
end
function PhotonStateful{Direction}(
mom::SFourMomentum
) where {Direction<:ParticleDirection}
return new{Direction,AbstractDefinitePolarization}(mom)
end
function PhotonStateful{Direction}(
mom::SFourMomentum, pol::AbstractDefinitePolarization
) where {Direction<:ParticleDirection}
return new{Direction,typeof(pol)}(mom)
end
function PhotonStateful{Dir,Pol}(ph::PhotonStateful) where {Dir,Pol}
return new{Dir,Pol}(ph.momentum)
end
end
PhotonStateful{Direction}(mom::SFourMomentum) where {Direction <: ParticleDirection} =
PhotonStateful{Direction, PolX}(mom)
PhotonStateful{Dir, Pol}(ph::PhotonStateful) where {Dir, Pol} = PhotonStateful{Dir, Pol}(ph.momentum)
"""
FermionStateful <: QEDParticle
A fermion of the [`QEDModel`](@ref) with its state.
"""
struct FermionStateful{Direction <: ParticleDirection, Spin <: AbstractDefiniteSpin} <: QEDParticle{Direction}
struct FermionStateful{Direction<:ParticleDirection,Spin<:AbstractDefiniteSpin} <:
QEDParticle{Direction}
momentum::SFourMomentum
# TODO: mass for electron/muon/tauon representation?
function FermionStateful{Direction,Spin}(
mom::SFourMomentum
) where {Direction<:ParticleDirection,Spin<:AbstractDefiniteSpin}
return new{Direction,Spin}(mom)
end
function FermionStateful{Direction}(
mom::SFourMomentum
) where {Direction<:ParticleDirection}
return new{Direction,AbstractDefiniteSpin}(mom)
end
function FermionStateful{Direction}(
mom::SFourMomentum, spin::AbstractDefiniteSpin
) where {Direction<:ParticleDirection}
return new{Direction,typeof(spin)}(mom)
end
function FermionStateful{Dir,Spin}(f::FermionStateful) where {Dir,Spin}
return new{Dir,Spin}(f.momentum)
end
end
FermionStateful{Direction}(mom::SFourMomentum) where {Direction <: ParticleDirection} =
FermionStateful{Direction, SpinUp}(mom)
FermionStateful{Dir, Spin}(f::FermionStateful) where {Dir, Spin} = FermionStateful{Dir, Spin}(f.momentum)
"""
AntiFermionStateful <: QEDParticle
An anti-fermion of the [`QEDModel`](@ref) with its state.
"""
struct AntiFermionStateful{Direction <: ParticleDirection, Spin <: AbstractDefiniteSpin} <: QEDParticle{Direction}
struct AntiFermionStateful{Direction<:ParticleDirection,Spin<:AbstractDefiniteSpin} <:
QEDParticle{Direction}
momentum::SFourMomentum
# TODO: mass for electron/muon/tauon representation?
function AntiFermionStateful{Direction,Spin}(
mom::SFourMomentum
) where {Direction<:ParticleDirection,Spin<:AbstractDefiniteSpin}
return new{Direction,Spin}(mom)
end
function AntiFermionStateful{Direction}(
mom::SFourMomentum
) where {Direction<:ParticleDirection}
return new{Direction,AbstractDefiniteSpin}(mom)
end
function AntiFermionStateful{Direction}(
mom::SFourMomentum, spin::AbstractDefiniteSpin
) where {Direction<:ParticleDirection}
return new{Direction,typeof(spin)}(mom)
end
function AntiFermionStateful{Dir,Spin}(f::AntiFermionStateful) where {Dir,Spin}
return new{Dir,Spin}(f.momentum)
end
end
AntiFermionStateful{Direction}(mom::SFourMomentum) where {Direction <: ParticleDirection} =
AntiFermionStateful{Direction, SpinUp}(mom)
AntiFermionStateful{Dir, Spin}(f::AntiFermionStateful) where {Dir, Spin} = AntiFermionStateful{Dir, Spin}(f.momentum)
"""
interaction_result(t1::Type{T1}, t2::Type{T2}) where {T1 <: QEDParticle, T2 <: QEDParticle}
For two given particle types that can interact, return the third.
"""
function interaction_result(t1::Type{T1}, t2::Type{T2}) where {T1 <: QEDParticle, T2 <: QEDParticle}
function interaction_result(
t1::Type{T1}, t2::Type{T2}
) where {T1<:QEDParticle,T2<:QEDParticle}
@assert false "Invalid interaction between particles of types $t1 and $t2"
end
interaction_result(
::Type{FermionStateful{Incoming, Spin1}},
::Type{FermionStateful{Outgoing, Spin2}},
) where {Spin1, Spin2} = PhotonStateful{Incoming, PolX}
interaction_result(
::Type{FermionStateful{Incoming, Spin1}},
::Type{AntiFermionStateful{Incoming, Spin2}},
) where {Spin1, Spin2} = PhotonStateful{Incoming, PolX}
interaction_result(::Type{FermionStateful{Incoming, Spin1}}, ::Type{<:PhotonStateful}) where {Spin1} =
FermionStateful{Outgoing, SpinUp}
function interaction_result(
::Type{FermionStateful{Incoming,Spin1}}, ::Type{FermionStateful{Outgoing,Spin2}}
) where {Spin1,Spin2}
return PhotonStateful{Incoming,PolX}
end
function interaction_result(
::Type{FermionStateful{Incoming,Spin1}}, ::Type{AntiFermionStateful{Incoming,Spin2}}
) where {Spin1,Spin2}
return PhotonStateful{Incoming,PolX}
end
function interaction_result(
::Type{FermionStateful{Incoming,Spin1}}, ::Type{<:PhotonStateful}
) where {Spin1}
return FermionStateful{Outgoing,SpinUp}
end
interaction_result(
::Type{FermionStateful{Outgoing, Spin1}},
::Type{FermionStateful{Incoming, Spin2}},
) where {Spin1, Spin2} = PhotonStateful{Incoming, PolX}
interaction_result(
::Type{FermionStateful{Outgoing, Spin1}},
::Type{AntiFermionStateful{Outgoing, Spin2}},
) where {Spin1, Spin2} = PhotonStateful{Incoming, PolX}
interaction_result(::Type{FermionStateful{Outgoing, Spin1}}, ::Type{<:PhotonStateful}) where {Spin1} =
FermionStateful{Incoming, SpinUp}
function interaction_result(
::Type{FermionStateful{Outgoing,Spin1}}, ::Type{FermionStateful{Incoming,Spin2}}
) where {Spin1,Spin2}
return PhotonStateful{Incoming,PolX}
end
function interaction_result(
::Type{FermionStateful{Outgoing,Spin1}}, ::Type{AntiFermionStateful{Outgoing,Spin2}}
) where {Spin1,Spin2}
return PhotonStateful{Incoming,PolX}
end
function interaction_result(
::Type{FermionStateful{Outgoing,Spin1}}, ::Type{<:PhotonStateful}
) where {Spin1}
return FermionStateful{Incoming,SpinUp}
end
# antifermion mirror
interaction_result(::Type{AntiFermionStateful{Incoming, Spin}}, t2::Type{<:QEDParticle}) where {Spin} =
interaction_result(FermionStateful{Outgoing, Spin}, t2)
interaction_result(::Type{AntiFermionStateful{Outgoing, Spin}}, t2::Type{<:QEDParticle}) where {Spin} =
interaction_result(FermionStateful{Incoming, Spin}, t2)
function interaction_result(
::Type{AntiFermionStateful{Incoming,Spin}}, t2::Type{<:QEDParticle}
) where {Spin}
return interaction_result(FermionStateful{Outgoing,Spin}, t2)
end
function interaction_result(
::Type{AntiFermionStateful{Outgoing,Spin}}, t2::Type{<:QEDParticle}
) where {Spin}
return interaction_result(FermionStateful{Incoming,Spin}, t2)
end
# photon commutativity
interaction_result(t1::Type{<:PhotonStateful}, t2::Type{<:QEDParticle}) = interaction_result(t2, t1)
function interaction_result(t1::Type{<:PhotonStateful}, t2::Type{<:QEDParticle})
return interaction_result(t2, t1)
end
# but prevent stack overflow
function interaction_result(t1::Type{<:PhotonStateful}, t2::Type{<:PhotonStateful})
@ -137,18 +211,34 @@ end
Return the type of the inverted direction. E.g.
"""
propagation_result(::Type{FermionStateful{Incoming, Spin}}) where {Spin <: AbstractDefiniteSpin} =
FermionStateful{Outgoing, Spin}
propagation_result(::Type{FermionStateful{Outgoing, Spin}}) where {Spin <: AbstractDefiniteSpin} =
FermionStateful{Incoming, Spin}
propagation_result(::Type{AntiFermionStateful{Incoming, Spin}}) where {Spin <: AbstractDefiniteSpin} =
AntiFermionStateful{Outgoing, Spin}
propagation_result(::Type{AntiFermionStateful{Outgoing, Spin}}) where {Spin <: AbstractDefiniteSpin} =
AntiFermionStateful{Incoming, Spin}
propagation_result(::Type{PhotonStateful{Incoming, Pol}}) where {Pol <: AbstractDefinitePolarization} =
PhotonStateful{Outgoing, Pol}
propagation_result(::Type{PhotonStateful{Outgoing, Pol}}) where {Pol <: AbstractDefinitePolarization} =
PhotonStateful{Incoming, Pol}
propagation_result(
::Type{FermionStateful{Incoming,Spin}}
) where {Spin<:AbstractDefiniteSpin} = FermionStateful{Outgoing,Spin}
function propagation_result(
::Type{FermionStateful{Outgoing,Spin}}
) where {Spin<:AbstractDefiniteSpin}
return FermionStateful{Incoming,Spin}
end
function propagation_result(
::Type{AntiFermionStateful{Incoming,Spin}}
) where {Spin<:AbstractDefiniteSpin}
return AntiFermionStateful{Outgoing,Spin}
end
function propagation_result(
::Type{AntiFermionStateful{Outgoing,Spin}}
) where {Spin<:AbstractDefiniteSpin}
return AntiFermionStateful{Incoming,Spin}
end
function propagation_result(
::Type{PhotonStateful{Incoming,Pol}}
) where {Pol<:AbstractDefinitePolarization}
return PhotonStateful{Outgoing,Pol}
end
function propagation_result(
::Type{PhotonStateful{Outgoing,Pol}}
) where {Pol<:AbstractDefinitePolarization}
return PhotonStateful{Incoming,Pol}
end
"""
types(::QEDModel)
@ -157,12 +247,12 @@ Return a Vector of the possible types of particle in the [`QEDModel`](@ref).
"""
function types(::QEDModel)
return [
PhotonStateful{Incoming, PolX},
PhotonStateful{Outgoing, PolX},
FermionStateful{Incoming, SpinUp},
FermionStateful{Outgoing, SpinUp},
AntiFermionStateful{Incoming, SpinUp},
AntiFermionStateful{Outgoing, SpinUp},
PhotonStateful{Incoming},
PhotonStateful{Outgoing},
FermionStateful{Incoming},
FermionStateful{Outgoing},
AntiFermionStateful{Incoming},
AntiFermionStateful{Outgoing},
]
end
@ -189,13 +279,13 @@ function String(::Type{<:AntiFermionStateful})
return "p"
end
function unique_name(::Type{PhotonStateful{Dir, Pol}}) where {Dir, Pol}
function unique_name(::Type{PhotonStateful{Dir,Pol}}) where {Dir,Pol}
return String(PhotonStateful) * String(Dir) * String(Pol)
end
function unique_name(::Type{FermionStateful{Dir, Spin}}) where {Dir, Spin}
function unique_name(::Type{FermionStateful{Dir,Spin}}) where {Dir,Spin}
return String(FermionStateful) * String(Dir) * String(Spin)
end
function unique_name(::Type{AntiFermionStateful{Dir, Spin}}) where {Dir, Spin}
function unique_name(::Type{AntiFermionStateful{Dir,Spin}}) where {Dir,Spin}
return String(AntiFermionStateful) * String(Dir) * String(Spin)
end
@ -203,27 +293,52 @@ end
@inline particle(::FermionStateful) = Electron()
@inline particle(::AntiFermionStateful) = Positron()
@inline particle(::Type{<:PhotonStateful}) = Photon()
@inline particle(::Type{<:FermionStateful}) = Electron()
@inline particle(::Type{<:AntiFermionStateful}) = Positron()
@inline momentum(p::PhotonStateful)::SFourMomentum = p.momentum
@inline momentum(p::FermionStateful)::SFourMomentum = p.momentum
@inline momentum(p::AntiFermionStateful)::SFourMomentum = p.momentum
@inline spin_or_pol(p::PhotonStateful{Dir, Pol}) where {Dir, Pol <: AbstractDefinitePolarization} = Pol()
@inline spin_or_pol(p::FermionStateful{Dir, Spin}) where {Dir, Spin <: AbstractDefiniteSpin} = Spin()
@inline spin_or_pol(p::AntiFermionStateful{Dir, Spin}) where {Dir, Spin <: AbstractDefiniteSpin} = Spin()
@inline spin_or_pol(
p::PhotonStateful{Dir,Pol}
) where {Dir,Pol<:AbstractDefinitePolarization} = Pol()
@inline spin_or_pol(p::FermionStateful{Dir,Spin}) where {Dir,Spin<:AbstractDefiniteSpin} =
Spin()
@inline spin_or_pol(
p::AntiFermionStateful{Dir,Spin}
) where {Dir,Spin<:AbstractDefiniteSpin} = Spin()
@inline direction(
::Type{P},
) where {P <: Union{FermionStateful{Incoming}, AntiFermionStateful{Incoming}, PhotonStateful{Incoming}}} = Incoming()
::Type{P}
) where {
P<:Union{
FermionStateful{Incoming},AntiFermionStateful{Incoming},PhotonStateful{Incoming}
},
} = Incoming()
@inline direction(
::Type{P},
) where {P <: Union{FermionStateful{Outgoing}, AntiFermionStateful{Outgoing}, PhotonStateful{Outgoing}}} = Outgoing()
::Type{P}
) where {
P<:Union{
FermionStateful{Outgoing},AntiFermionStateful{Outgoing},PhotonStateful{Outgoing}
},
} = Outgoing()
@inline direction(
::P,
) where {P <: Union{FermionStateful{Incoming}, AntiFermionStateful{Incoming}, PhotonStateful{Incoming}}} = Incoming()
::P
) where {
P<:Union{
FermionStateful{Incoming},AntiFermionStateful{Incoming},PhotonStateful{Incoming}
},
} = Incoming()
@inline direction(
::P,
) where {P <: Union{FermionStateful{Outgoing}, AntiFermionStateful{Outgoing}, PhotonStateful{Outgoing}}} = Outgoing()
::P
) where {
P<:Union{
FermionStateful{Outgoing},AntiFermionStateful{Outgoing},PhotonStateful{Outgoing}
},
} = Outgoing()
@inline isincoming(::QEDParticle{Incoming}) = true
@inline isincoming(::QEDParticle{Outgoing}) = false
@ -239,13 +354,12 @@ end
@inline mass(::Type{<:AntiFermionStateful}) = 1.0
@inline mass(::Type{<:PhotonStateful}) = 0.0
@inline invert_momentum(p::FermionStateful{Dir, Spin}) where {Dir, Spin} =
FermionStateful{Dir, Spin}(-p.momentum, p.spin)
@inline invert_momentum(p::AntiFermionStateful{Dir, Spin}) where {Dir, Spin} =
AntiFermionStateful{Dir, Spin}(-p.momentum, p.spin)
@inline invert_momentum(k::PhotonStateful{Dir, Spin}) where {Dir, Spin} =
PhotonStateful{Dir, Spin}(-k.momentum, k.polarization)
@inline invert_momentum(p::FermionStateful{Dir,Spin}) where {Dir,Spin} =
FermionStateful{Dir,Spin}(-p.momentum, p.spin)
@inline invert_momentum(p::AntiFermionStateful{Dir,Spin}) where {Dir,Spin} =
AntiFermionStateful{Dir,Spin}(-p.momentum, p.spin)
@inline invert_momentum(k::PhotonStateful{Dir,Spin}) where {Dir,Spin} =
PhotonStateful{Dir,Spin}(-k.momentum, k.polarization)
"""
caninteract(T1::Type{<:QEDParticle}, T2::Type{<:QEDParticle})
@ -276,17 +390,17 @@ end
function type_index_from_name(::QEDModel, name::String)
if startswith(name, "ki")
return (PhotonStateful{Incoming, PolX}, parse(Int, name[3:end]))
return (PhotonStateful{Incoming,PolX}, parse(Int, name[3:end]))
elseif startswith(name, "ko")
return (PhotonStateful{Outgoing, PolX}, parse(Int, name[3:end]))
return (PhotonStateful{Outgoing,PolX}, parse(Int, name[3:end]))
elseif startswith(name, "ei")
return (FermionStateful{Incoming, SpinUp}, parse(Int, name[3:end]))
return (FermionStateful{Incoming,SpinUp}, parse(Int, name[3:end]))
elseif startswith(name, "eo")
return (FermionStateful{Outgoing, SpinUp}, parse(Int, name[3:end]))
return (FermionStateful{Outgoing,SpinUp}, parse(Int, name[3:end]))
elseif startswith(name, "pi")
return (AntiFermionStateful{Incoming, SpinUp}, parse(Int, name[3:end]))
return (AntiFermionStateful{Incoming,SpinUp}, parse(Int, name[3:end]))
elseif startswith(name, "po")
return (AntiFermionStateful{Outgoing, SpinUp}, parse(Int, name[3:end]))
return (AntiFermionStateful{Outgoing,SpinUp}, parse(Int, name[3:end]))
else
throw("Invalid name for a particle in the QED model")
end
@ -314,7 +428,7 @@ Return the factor of a vertex in a QED feynman diagram.
end
@inline function QED_inner_edge(p::QEDParticle)
return propagator(particle(p), p.momentum)
return QEDbase.propagator(particle(p), p.momentum)
end
"""
@ -323,15 +437,22 @@ end
Calculate and return a new particle from two given interacting ones at a vertex.
"""
function QED_conserve_momentum(
p1::P1,
p2::P2,
p1::P1, p2::P2
) where {
Dir1 <: ParticleDirection,
Dir2 <: ParticleDirection,
SpinPol1 <: AbstractSpinOrPolarization,
SpinPol2 <: AbstractSpinOrPolarization,
P1 <: Union{FermionStateful{Dir1, SpinPol1}, AntiFermionStateful{Dir1, SpinPol1}, PhotonStateful{Dir1, SpinPol1}},
P2 <: Union{FermionStateful{Dir2, SpinPol2}, AntiFermionStateful{Dir2, SpinPol2}, PhotonStateful{Dir2, SpinPol2}},
Dir1<:ParticleDirection,
Dir2<:ParticleDirection,
SpinPol1<:AbstractSpinOrPolarization,
SpinPol2<:AbstractSpinOrPolarization,
P1<:Union{
FermionStateful{Dir1,SpinPol1},
AntiFermionStateful{Dir1,SpinPol1},
PhotonStateful{Dir1,SpinPol1},
},
P2<:Union{
FermionStateful{Dir2,SpinPol2},
AntiFermionStateful{Dir2,SpinPol2},
PhotonStateful{Dir2,SpinPol2},
},
}
P3 = interaction_result(P1, P2)
p1_mom = p1.momentum
@ -357,14 +478,14 @@ Input for a QED Process. Contains the [`QEDProcessDescription`](@ref) of the pro
See also: [`gen_process_input`](@ref)
"""
struct QEDProcessInput{N1, N2, N3, N4, N5, N6} <: AbstractProcessInput
struct QEDProcessInput{N1,N2,N3,N4,N5,N6,S1,S2,S3,S4,P1,P2} <: AbstractProcessInput
process::QEDProcessDescription
inFerms::SVector{N1, FermionStateful{Incoming, SpinUp}}
outFerms::SVector{N2, FermionStateful{Outgoing, SpinUp}}
inAntiferms::SVector{N3, AntiFermionStateful{Incoming, SpinUp}}
outAntiferms::SVector{N4, AntiFermionStateful{Outgoing, SpinUp}}
inPhotons::SVector{N5, PhotonStateful{Incoming, PolX}}
outPhotons::SVector{N6, PhotonStateful{Outgoing, PolX}}
inFerms::SVector{N1,FermionStateful{Incoming,S1}}
outFerms::SVector{N2,FermionStateful{Outgoing,S2}}
inAntiferms::SVector{N3,AntiFermionStateful{Incoming,S3}}
outAntiferms::SVector{N4,AntiFermionStateful{Outgoing,S4}}
inPhotons::SVector{N5,PhotonStateful{Incoming,P1}}
outPhotons::SVector{N6,PhotonStateful{Outgoing,P2}}
end
"""
@ -379,8 +500,9 @@ function copy(process::QEDProcessDescription)
return QEDProcessDescription(copy(process.inParticles), copy(process.outParticles))
end
==(p1::QEDProcessDescription, p2::QEDProcessDescription) =
p1.inParticles == p2.inParticles && p1.outParticles == p2.outParticles
function ==(p1::QEDProcessDescription, p2::QEDProcessDescription)
return p1.inParticles == p2.inParticles && p1.outParticles == p2.outParticles
end
function in_particles(process::QEDProcessDescription)
return process.inParticles
@ -390,7 +512,9 @@ function out_particles(process::QEDProcessDescription)
return process.outParticles
end
function get_particle(input::QEDProcessInput, t::Type{Particle}, n::Int)::Particle where {Particle}
function get_particle(
input::QEDProcessInput, t::Type{Particle}, n::Int
)::Particle where {Particle}
if (t <: FermionStateful{Incoming})
return input.inFerms[n]
elseif (t <: FermionStateful{Outgoing})

View File

@ -18,21 +18,21 @@ def_momentum = SFourMomentum(1.0, 0.0, 0.0, 0.0)
RNG = Random.MersenneTwister(0)
testparticleTypes = [
PhotonStateful{Incoming, PolX},
PhotonStateful{Outgoing, PolX},
FermionStateful{Incoming, SpinUp},
FermionStateful{Outgoing, SpinUp},
AntiFermionStateful{Incoming, SpinUp},
AntiFermionStateful{Outgoing, SpinUp},
PhotonStateful{Incoming,PolX},
PhotonStateful{Outgoing,PolX},
FermionStateful{Incoming,SpinUp},
FermionStateful{Outgoing,SpinUp},
AntiFermionStateful{Incoming,SpinUp},
AntiFermionStateful{Outgoing,SpinUp},
]
testparticleTypesPropagated = [
PhotonStateful{Outgoing, PolX},
PhotonStateful{Incoming, PolX},
FermionStateful{Outgoing, SpinUp},
FermionStateful{Incoming, SpinUp},
AntiFermionStateful{Outgoing, SpinUp},
AntiFermionStateful{Incoming, SpinUp},
PhotonStateful{Outgoing,PolX},
PhotonStateful{Incoming,PolX},
FermionStateful{Outgoing,SpinUp},
FermionStateful{Incoming,SpinUp},
AntiFermionStateful{Outgoing,SpinUp},
AntiFermionStateful{Incoming,SpinUp},
]
function compton_groundtruth(input::QEDProcessInput)
@ -57,8 +57,8 @@ function compton_groundtruth(input::QEDProcessInput)
virt2_mom = p1.momentum + k1.momentum
@test isapprox(p2.momentum + k2.momentum, virt2_mom)
s_p2_k1 = propagator(Electron(), virt1_mom)
s_p1_k1 = propagator(Electron(), virt2_mom)
s_p2_k1 = QEDbase.propagator(Electron(), virt1_mom)
s_p1_k1 = QEDbase.propagator(Electron(), virt2_mom)
diagram1 = u_p2 * (eps_1 * QED_vertex()) * s_p2_k1 * (eps_2 * QED_vertex()) * u_p1
diagram2 = u_p2 * (eps_2 * QED_vertex()) * s_p1_k1 * (eps_1 * QED_vertex()) * u_p1
@ -66,7 +66,6 @@ function compton_groundtruth(input::QEDProcessInput)
return diagram1 + diagram2
end
@testset "Interaction Result" begin
import MetagraphOptimization.QED_conserve_momentum
@ -88,7 +87,11 @@ end
@test issame(typeof(resultParticle), interaction_result(p1, p2))
totalMom = zero(SFourMomentum)
for (p, mom) in [(p1, testParticle1.momentum), (p2, testParticle2.momentum), (p3, resultParticle.momentum)]
for (p, mom) in [
(p1, testParticle1.momentum),
(p2, testParticle2.momentum),
(p3, resultParticle.momentum),
]
if (typeof(direction(p)) <: Incoming)
totalMom += mom
else
@ -96,7 +99,7 @@ end
end
end
@test isapprox(totalMom, zero(SFourMomentum); atol = sqrt(eps()))
@test isapprox(totalMom, zero(SFourMomentum); atol=sqrt(eps()))
end
end
@ -113,41 +116,63 @@ end
@test parse_process("ke->ke", QEDModel()) == parse_process("ek->ek", QEDModel())
@test parse_process("ke->ke", QEDModel()) == parse_process("ke->ek", QEDModel())
@test parse_process("kkke->eep", QEDModel()) == parse_process("kkek->epe", QEDModel())
@test parse_process("kkke->eep", QEDModel()) ==
parse_process("kkek->epe", QEDModel())
end
@testset "Known processes" begin
compton_process = QEDProcessDescription(
Dict{Type, Int}(PhotonStateful{Incoming, PolX} => 1, FermionStateful{Incoming, SpinUp} => 1),
Dict{Type, Int}(PhotonStateful{Outgoing, PolX} => 1, FermionStateful{Outgoing, SpinUp} => 1),
Dict{Type,Int}(
PhotonStateful{Incoming,PolX} => 1, FermionStateful{Incoming,SpinUp} => 1
),
Dict{Type,Int}(
PhotonStateful{Outgoing,PolX} => 1, FermionStateful{Outgoing,SpinUp} => 1
),
)
@test parse_process("ke->ke", QEDModel()) == compton_process
positron_compton_process = QEDProcessDescription(
Dict{Type, Int}(PhotonStateful{Incoming, PolX} => 1, AntiFermionStateful{Incoming, SpinUp} => 1),
Dict{Type, Int}(PhotonStateful{Outgoing, PolX} => 1, AntiFermionStateful{Outgoing, SpinUp} => 1),
Dict{Type,Int}(
PhotonStateful{Incoming,PolX} => 1,
AntiFermionStateful{Incoming,SpinUp} => 1,
),
Dict{Type,Int}(
PhotonStateful{Outgoing,PolX} => 1,
AntiFermionStateful{Outgoing,SpinUp} => 1,
),
)
@test parse_process("kp->kp", QEDModel()) == positron_compton_process
trident_process = QEDProcessDescription(
Dict{Type, Int}(PhotonStateful{Incoming, PolX} => 1, FermionStateful{Incoming, SpinUp} => 1),
Dict{Type, Int}(FermionStateful{Outgoing, SpinUp} => 2, AntiFermionStateful{Outgoing, SpinUp} => 1),
Dict{Type,Int}(
PhotonStateful{Incoming,PolX} => 1, FermionStateful{Incoming,SpinUp} => 1
),
Dict{Type,Int}(
FermionStateful{Outgoing,SpinUp} => 2,
AntiFermionStateful{Outgoing,SpinUp} => 1,
),
)
@test parse_process("ke->eep", QEDModel()) == trident_process
pair_production_process = QEDProcessDescription(
Dict{Type, Int}(PhotonStateful{Incoming, PolX} => 2),
Dict{Type, Int}(FermionStateful{Outgoing, SpinUp} => 1, AntiFermionStateful{Outgoing, SpinUp} => 1),
Dict{Type,Int}(PhotonStateful{Incoming,PolX} => 2),
Dict{Type,Int}(
FermionStateful{Outgoing,SpinUp} => 1,
AntiFermionStateful{Outgoing,SpinUp} => 1,
),
)
@test parse_process("kk->pe", QEDModel()) == pair_production_process
pair_annihilation_process = QEDProcessDescription(
Dict{Type, Int}(FermionStateful{Incoming, SpinUp} => 1, AntiFermionStateful{Incoming, SpinUp} => 1),
Dict{Type, Int}(PhotonStateful{Outgoing, PolX} => 2),
Dict{Type,Int}(
FermionStateful{Incoming,SpinUp} => 1,
AntiFermionStateful{Incoming,SpinUp} => 1,
),
Dict{Type,Int}(PhotonStateful{Outgoing,PolX} => 2),
)
@test parse_process("pe->kk", QEDModel()) == pair_annihilation_process
@ -161,12 +186,18 @@ end
for i in 1:100
input = gen_process_input(process)
@test length(input.inFerms) == get(process.inParticles, FermionStateful{Incoming, SpinUp}, 0)
@test length(input.inAntiferms) == get(process.inParticles, AntiFermionStateful{Incoming, SpinUp}, 0)
@test length(input.inPhotons) == get(process.inParticles, PhotonStateful{Incoming, PolX}, 0)
@test length(input.outFerms) == get(process.outParticles, FermionStateful{Outgoing, SpinUp}, 0)
@test length(input.outAntiferms) == get(process.outParticles, AntiFermionStateful{Outgoing, SpinUp}, 0)
@test length(input.outPhotons) == get(process.outParticles, PhotonStateful{Outgoing, PolX}, 0)
@test length(input.inFerms) ==
get(process.inParticles, FermionStateful{Incoming,SpinUp}, 0)
@test length(input.inAntiferms) ==
get(process.inParticles, AntiFermionStateful{Incoming,SpinUp}, 0)
@test length(input.inPhotons) ==
get(process.inParticles, PhotonStateful{Incoming,PolX}, 0)
@test length(input.outFerms) ==
get(process.outParticles, FermionStateful{Outgoing,SpinUp}, 0)
@test length(input.outAntiferms) ==
get(process.outParticles, AntiFermionStateful{Outgoing,SpinUp}, 0)
@test length(input.outPhotons) ==
get(process.outParticles, PhotonStateful{Outgoing,PolX}, 0)
@test isapprox(
sum([
@ -179,7 +210,7 @@ end
getfield.(input.outAntiferms, :momentum)...,
getfield.(input.outPhotons, :momentum)...,
]);
atol = sqrt(eps()),
atol=sqrt(eps()),
)
end
end
@ -211,97 +242,97 @@ end
graph = DAG()
# s to output (exit node)
d_exit = insert_node!(graph, make_node(DataTask(16)), track = false)
d_exit = insert_node!(graph, make_node(DataTask(16)); track=false)
sum_node = insert_node!(graph, make_node(ComputeTaskQED_Sum(2)), track = false)
sum_node = insert_node!(graph, make_node(ComputeTaskQED_Sum(2)); track=false)
d_s0_sum = insert_node!(graph, make_node(DataTask(16)), track = false)
d_s1_sum = insert_node!(graph, make_node(DataTask(16)), track = false)
d_s0_sum = insert_node!(graph, make_node(DataTask(16)); track=false)
d_s1_sum = insert_node!(graph, make_node(DataTask(16)); track=false)
# final s compute
s0 = insert_node!(graph, make_node(ComputeTaskQED_S2()), track = false)
s1 = insert_node!(graph, make_node(ComputeTaskQED_S2()), track = false)
s0 = insert_node!(graph, make_node(ComputeTaskQED_S2()); track=false)
s1 = insert_node!(graph, make_node(ComputeTaskQED_S2()); track=false)
# data from v0 and v1 to s0
d_v0_s0 = insert_node!(graph, make_node(DataTask(96)), track = false)
d_v1_s0 = insert_node!(graph, make_node(DataTask(96)), track = false)
d_v2_s1 = insert_node!(graph, make_node(DataTask(96)), track = false)
d_v3_s1 = insert_node!(graph, make_node(DataTask(96)), track = false)
d_v0_s0 = insert_node!(graph, make_node(DataTask(96)); track=false)
d_v1_s0 = insert_node!(graph, make_node(DataTask(96)); track=false)
d_v2_s1 = insert_node!(graph, make_node(DataTask(96)); track=false)
d_v3_s1 = insert_node!(graph, make_node(DataTask(96)); track=false)
# v0 and v1 compute
v0 = insert_node!(graph, make_node(ComputeTaskQED_V()), track = false)
v1 = insert_node!(graph, make_node(ComputeTaskQED_V()), track = false)
v2 = insert_node!(graph, make_node(ComputeTaskQED_V()), track = false)
v3 = insert_node!(graph, make_node(ComputeTaskQED_V()), track = false)
v0 = insert_node!(graph, make_node(ComputeTaskQED_V()); track=false)
v1 = insert_node!(graph, make_node(ComputeTaskQED_V()); track=false)
v2 = insert_node!(graph, make_node(ComputeTaskQED_V()); track=false)
v3 = insert_node!(graph, make_node(ComputeTaskQED_V()); track=false)
# data from uPhIn, uPhOut, uElIn, uElOut to v0 and v1
d_uPhIn_v0 = insert_node!(graph, make_node(DataTask(96)), track = false)
d_uElIn_v0 = insert_node!(graph, make_node(DataTask(96)), track = false)
d_uPhOut_v1 = insert_node!(graph, make_node(DataTask(96)), track = false)
d_uElOut_v1 = insert_node!(graph, make_node(DataTask(96)), track = false)
d_uPhIn_v0 = insert_node!(graph, make_node(DataTask(96)); track=false)
d_uElIn_v0 = insert_node!(graph, make_node(DataTask(96)); track=false)
d_uPhOut_v1 = insert_node!(graph, make_node(DataTask(96)); track=false)
d_uElOut_v1 = insert_node!(graph, make_node(DataTask(96)); track=false)
# data from uPhIn, uPhOut, uElIn, uElOut to v2 and v3
d_uPhOut_v2 = insert_node!(graph, make_node(DataTask(96)), track = false)
d_uElIn_v2 = insert_node!(graph, make_node(DataTask(96)), track = false)
d_uPhIn_v3 = insert_node!(graph, make_node(DataTask(96)), track = false)
d_uElOut_v3 = insert_node!(graph, make_node(DataTask(96)), track = false)
d_uPhOut_v2 = insert_node!(graph, make_node(DataTask(96)); track=false)
d_uElIn_v2 = insert_node!(graph, make_node(DataTask(96)); track=false)
d_uPhIn_v3 = insert_node!(graph, make_node(DataTask(96)); track=false)
d_uElOut_v3 = insert_node!(graph, make_node(DataTask(96)); track=false)
# uPhIn, uPhOut, uElIn and uElOut computes
uPhIn = insert_node!(graph, make_node(ComputeTaskQED_U()), track = false)
uPhOut = insert_node!(graph, make_node(ComputeTaskQED_U()), track = false)
uElIn = insert_node!(graph, make_node(ComputeTaskQED_U()), track = false)
uElOut = insert_node!(graph, make_node(ComputeTaskQED_U()), track = false)
uPhIn = insert_node!(graph, make_node(ComputeTaskQED_U()); track=false)
uPhOut = insert_node!(graph, make_node(ComputeTaskQED_U()); track=false)
uElIn = insert_node!(graph, make_node(ComputeTaskQED_U()); track=false)
uElOut = insert_node!(graph, make_node(ComputeTaskQED_U()); track=false)
# data into U
d_uPhIn = insert_node!(graph, make_node(DataTask(16), "ki1"), track = false)
d_uPhOut = insert_node!(graph, make_node(DataTask(16), "ko1"), track = false)
d_uElIn = insert_node!(graph, make_node(DataTask(16), "ei1"), track = false)
d_uElOut = insert_node!(graph, make_node(DataTask(16), "eo1"), track = false)
d_uPhIn = insert_node!(graph, make_node(DataTask(16), "ki1"); track=false)
d_uPhOut = insert_node!(graph, make_node(DataTask(16), "ko1"); track=false)
d_uElIn = insert_node!(graph, make_node(DataTask(16), "ei1"); track=false)
d_uElOut = insert_node!(graph, make_node(DataTask(16), "eo1"); track=false)
# now for all the edges
insert_edge!(graph, d_uPhIn, uPhIn, track = false)
insert_edge!(graph, d_uPhOut, uPhOut, track = false)
insert_edge!(graph, d_uElIn, uElIn, track = false)
insert_edge!(graph, d_uElOut, uElOut, track = false)
insert_edge!(graph, d_uPhIn, uPhIn; track=false)
insert_edge!(graph, d_uPhOut, uPhOut; track=false)
insert_edge!(graph, d_uElIn, uElIn; track=false)
insert_edge!(graph, d_uElOut, uElOut; track=false)
insert_edge!(graph, uPhIn, d_uPhIn_v0, track = false)
insert_edge!(graph, uPhOut, d_uPhOut_v1, track = false)
insert_edge!(graph, uElIn, d_uElIn_v0, track = false)
insert_edge!(graph, uElOut, d_uElOut_v1, track = false)
insert_edge!(graph, uPhIn, d_uPhIn_v0; track=false)
insert_edge!(graph, uPhOut, d_uPhOut_v1; track=false)
insert_edge!(graph, uElIn, d_uElIn_v0; track=false)
insert_edge!(graph, uElOut, d_uElOut_v1; track=false)
insert_edge!(graph, uPhIn, d_uPhIn_v3, track = false)
insert_edge!(graph, uPhOut, d_uPhOut_v2, track = false)
insert_edge!(graph, uElIn, d_uElIn_v2, track = false)
insert_edge!(graph, uElOut, d_uElOut_v3, track = false)
insert_edge!(graph, uPhIn, d_uPhIn_v3; track=false)
insert_edge!(graph, uPhOut, d_uPhOut_v2; track=false)
insert_edge!(graph, uElIn, d_uElIn_v2; track=false)
insert_edge!(graph, uElOut, d_uElOut_v3; track=false)
insert_edge!(graph, d_uPhIn_v0, v0, track = false)
insert_edge!(graph, d_uPhOut_v1, v1, track = false)
insert_edge!(graph, d_uElIn_v0, v0, track = false)
insert_edge!(graph, d_uElOut_v1, v1, track = false)
insert_edge!(graph, d_uPhIn_v0, v0; track=false)
insert_edge!(graph, d_uPhOut_v1, v1; track=false)
insert_edge!(graph, d_uElIn_v0, v0; track=false)
insert_edge!(graph, d_uElOut_v1, v1; track=false)
insert_edge!(graph, d_uPhIn_v3, v3, track = false)
insert_edge!(graph, d_uPhOut_v2, v2, track = false)
insert_edge!(graph, d_uElIn_v2, v2, track = false)
insert_edge!(graph, d_uElOut_v3, v3, track = false)
insert_edge!(graph, d_uPhIn_v3, v3; track=false)
insert_edge!(graph, d_uPhOut_v2, v2; track=false)
insert_edge!(graph, d_uElIn_v2, v2; track=false)
insert_edge!(graph, d_uElOut_v3, v3; track=false)
insert_edge!(graph, v0, d_v0_s0, track = false)
insert_edge!(graph, v1, d_v1_s0, track = false)
insert_edge!(graph, v2, d_v2_s1, track = false)
insert_edge!(graph, v3, d_v3_s1, track = false)
insert_edge!(graph, v0, d_v0_s0; track=false)
insert_edge!(graph, v1, d_v1_s0; track=false)
insert_edge!(graph, v2, d_v2_s1; track=false)
insert_edge!(graph, v3, d_v3_s1; track=false)
insert_edge!(graph, d_v0_s0, s0, track = false)
insert_edge!(graph, d_v1_s0, s0, track = false)
insert_edge!(graph, d_v0_s0, s0; track=false)
insert_edge!(graph, d_v1_s0, s0; track=false)
insert_edge!(graph, d_v2_s1, s1, track = false)
insert_edge!(graph, d_v3_s1, s1, track = false)
insert_edge!(graph, d_v2_s1, s1; track=false)
insert_edge!(graph, d_v3_s1, s1; track=false)
insert_edge!(graph, s0, d_s0_sum, track = false)
insert_edge!(graph, s1, d_s1_sum, track = false)
insert_edge!(graph, s0, d_s0_sum; track=false)
insert_edge!(graph, s1, d_s1_sum; track=false)
insert_edge!(graph, d_s0_sum, sum_node, track = false)
insert_edge!(graph, d_s1_sum, sum_node, track = false)
insert_edge!(graph, d_s0_sum, sum_node; track=false)
insert_edge!(graph, d_s1_sum, sum_node; track=false)
insert_edge!(graph, sum_node, d_exit, track = false)
insert_edge!(graph, sum_node, d_exit; track=false)
input = [gen_process_input(process) for _ in 1:1000]
@ -314,9 +345,12 @@ end
@test isapprox(compton_function.(input), compton_groundtruth.(input))
end
@testset "Equal results after optimization" for optimizer in
[ReductionOptimizer(), RandomWalkOptimizer(MersenneTwister(0))]
@testset "Process $proc_str" for proc_str in ["ke->ke", "kp->kp", "kk->ep", "ep->kk", "ke->kke", "ke->kkke"]
@testset "Equal results after optimization" for optimizer in [
ReductionOptimizer(), RandomWalkOptimizer(MersenneTwister(0))
]
@testset "Process $proc_str" for proc_str in [
"ke->ke", "kp->kp", "kk->ep", "ep->kk", "ke->kke", "ke->kkke"
]
model = QEDModel()
process = parse_process(proc_str, model)
machine = Machine(