From 901944bd8b5bd617859aa59f89b61fb4069a03e5 Mon Sep 17 00:00:00 2001 From: Rubydragon Date: Thu, 27 Jun 2024 13:21:53 +0200 Subject: [PATCH] WIP --- Project.toml | 1 + examples/congruent_in_ph.jl | 123 ++++++++++++ src/MetagraphOptimization.jl | 3 +- src/models/qed/compute.jl | 37 ++-- src/models/qed/create.jl | 162 +++++++++++----- src/models/qed/diagrams.jl | 127 +++++++------ src/models/qed/particle.jl | 358 +++++++++++++++++++++++------------ test/unit_tests_qedmodel.jl | 236 +++++++++++++---------- 8 files changed, 710 insertions(+), 337 deletions(-) create mode 100644 examples/congruent_in_ph.jl diff --git a/Project.toml b/Project.toml index 63f9d2e..32c8b7b 100644 --- a/Project.toml +++ b/Project.toml @@ -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" diff --git a/examples/congruent_in_ph.jl b/examples/congruent_in_ph.jl new file mode 100644 index 0000000..a0639bf --- /dev/null +++ b/examples/congruent_in_ph.jl @@ -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 diff --git a/src/MetagraphOptimization.jl b/src/MetagraphOptimization.jl index 608b515..6240955 100644 --- a/src/MetagraphOptimization.jl +++ b/src/MetagraphOptimization.jl @@ -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") diff --git a/src/models/qed/compute.jl b/src/models/qed/compute.jl index b30994f..ad4a5a6 100644 --- a/src/models/qed/compute.jl +++ b/src/models/qed/compute.jl @@ -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 diff --git a/src/models/qed/create.jl b/src/models/qed/create.jl index 85ee3aa..b853d3f 100644 --- a/src/models/qed/create.jl +++ b/src/models/qed/create.jl @@ -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 diff --git a/src/models/qed/diagrams.jl b/src/models/qed/diagrams.jl index 9477eb6..6c882de 100644 --- a/src/models/qed/diagrams.jl +++ b/src/models/qed/diagrams.jl @@ -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 diff --git a/src/models/qed/particle.jl b/src/models/qed/particle.jl index d51fedd..8e05107 100644 --- a/src/models/qed/particle.jl +++ b/src/models/qed/particle.jl @@ -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}) diff --git a/test/unit_tests_qedmodel.jl b/test/unit_tests_qedmodel.jl index 8126605..8df14d5 100644 --- a/test/unit_tests_qedmodel.jl +++ b/test/unit_tests_qedmodel.jl @@ -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(