From 29071a9cb23ee0945d3b6739d840af81f9dfcd14 Mon Sep 17 00:00:00 2001 From: Anton Reinhard Date: Wed, 19 Jun 2024 16:47:29 +0200 Subject: [PATCH] Remove process from diagram definition --- src/FeynmanDiagramGenerator.jl | 6 +++++ src/QEDprocesses_patch.jl | 7 ++++++ src/diagrams/diagrams.jl | 46 ++++++++++++++++++++-------------- src/generic_process_def.jl | 34 +++++++++++++++++++++++++ 4 files changed, 74 insertions(+), 19 deletions(-) create mode 100644 src/QEDprocesses_patch.jl create mode 100644 src/generic_process_def.jl diff --git a/src/FeynmanDiagramGenerator.jl b/src/FeynmanDiagramGenerator.jl index 450c5b1..a460378 100644 --- a/src/FeynmanDiagramGenerator.jl +++ b/src/FeynmanDiagramGenerator.jl @@ -3,8 +3,12 @@ module FeynmanDiagramGenerator using QEDbase using QEDprocesses +include("QEDprocesses_patch.jl") + import Base.== +export GenericQEDProcess, isphysical + export AbstractTreeLevelFeynmanDiagram, FeynmanVertex, FeynmanDiagram export external_particles, virtual_particles, process, vertices @@ -22,6 +26,8 @@ include("trees/print.jl") include("qft/qft.jl") +include("generic_process_def.jl") + include("diagrams/interface.jl") include("diagrams/vertex.jl") include("diagrams/diagrams.jl") diff --git a/src/QEDprocesses_patch.jl b/src/QEDprocesses_patch.jl new file mode 100644 index 0000000..ba08085 --- /dev/null +++ b/src/QEDprocesses_patch.jl @@ -0,0 +1,7 @@ +# patch QEDprocesses +# see issue https://github.com/QEDjl-project/QEDprocesses.jl/issues/77 +@inline function QEDprocesses.number_particles( + proc_def::AbstractProcessDefinition, dir::DIR, ::PT +) where {DIR<:ParticleDirection,PT<:AbstractParticleType} + return count(x -> x isa PT, particles(proc_def, dir)) +end diff --git a/src/diagrams/diagrams.jl b/src/diagrams/diagrams.jl index be79384..159ce44 100644 --- a/src/diagrams/diagrams.jl +++ b/src/diagrams/diagrams.jl @@ -43,10 +43,8 @@ end # Feynman Diagram, tree-level, QED # -struct FeynmanDiagram{PROC,N,E,U,T,M} <: AbstractTreeLevelFeynmanDiagram where {P<:AbstractProcessDefinition,N,E,U,T,M} - # E, U, T, and M can be inferred from the PROC, but not necessarily in 0 runtime - proc::PROC - +struct FeynmanDiagram{N,E,U,T,M} <: AbstractTreeLevelFeynmanDiagram where {N,E,U,T,M} + # TODO: flatten into one list diagram_structure::NTuple{N,Vector{Int}} electron_permutation::NTuple{E,Int} @@ -54,7 +52,6 @@ struct FeynmanDiagram{PROC,N,E,U,T,M} <: AbstractTreeLevelFeynmanDiagram where { tauon_permutation::NTuple{T,Int} function FeynmanDiagram( - proc::PROC, structure::Vector{Vector{Int}}, elec_perm::Vector{Int}, muon_perm::Vector{Int}, @@ -63,18 +60,16 @@ struct FeynmanDiagram{PROC,N,E,U,T,M} <: AbstractTreeLevelFeynmanDiagram where { ::Val{U}, ::Val{T}, ::Val{M} - ) where {PROC<:AbstractProcessDefinition,E,U,T,M} + ) where {E,U,T,M} @assert E == length(elec_perm) @assert U == length(muon_perm) @assert T == length(tauon_perm) N = E + U + T - return new{PROC,N,E,U,T,M}(proc, NTuple{N,Vector{Int}}(structure), NTuple{E,Int}(elec_perm), NTuple{U,Int}(muon_perm), NTuple{T,Int}(tauon_perm)) + return new{N,E,U,T,M}(NTuple{N,Vector{Int}}(structure), NTuple{E,Int}(elec_perm), NTuple{U,Int}(muon_perm), NTuple{T,Int}(tauon_perm)) end end -process(diagram::FeynmanDiagram) = d.proc - function virtual_particles(diagram::FeynmanDiagram) return NTuple{N,Tuple{QEDbase.AbstractParticleType,BitArray}}() @@ -159,8 +154,7 @@ function _external_photon_multiplicity(fermion_structure::Vector{Vector{Int}}, n return res end -mutable struct FeynmanDiagramIterator{PROC<:AbstractProcessDefinition,E,U,T,M} - process::PROC +mutable struct FeynmanDiagramIterator{E,U,T,M} e::Val{E} # number of electron lines (indices 1 - e) e_perms::Vector{Vector{Int}} # list of all the possible permutations of the electrons e_index::Int @@ -175,14 +169,14 @@ mutable struct FeynmanDiagramIterator{PROC<:AbstractProcessDefinition,E,U,T,M} photon_structure::Vector{Vector{Int}} # current structure that's being permuted end -function Base.length(it::FeynmanDiagramIterator{<:AbstractProcessDefinition,E,U,T,M}) where {E,U,T,M} +function Base.length(it::FeynmanDiagramIterator{E,U,T,M}) where {E,U,T,M} N = E + U + T return factorial(M + 3 * N - 3, 2 * N - 1) * factorial(E) * factorial(U) * factorial(T) end function Base.iterate(iter::FeynmanDiagramIterator) return ( - FeynmanDiagram(iter.process, iter.photon_structure, iter.e_perms[iter.e_index], iter.u_perms[iter.u_index], iter.t_perms[iter.t_index], iter.e, iter.u, iter.t, iter.m), + FeynmanDiagram(iter.photon_structure, iter.e_perms[iter.e_index], iter.u_perms[iter.u_index], iter.t_perms[iter.t_index], iter.e, iter.u, iter.t, iter.m), nothing ) end @@ -210,15 +204,29 @@ function Base.iterate(iter::FeynmanDiagramIterator, ::Nothing) end return ( - FeynmanDiagram(iter.process, iter.photon_structure, iter.e_perms[iter.e_index], iter.u_perms[iter.u_index], iter.t_perms[iter.t_index], iter.e, iter.u, iter.t, iter.m), + FeynmanDiagram(iter.photon_structure, iter.e_perms[iter.e_index], iter.u_perms[iter.u_index], iter.t_perms[iter.t_index], iter.e, iter.u, iter.t, iter.m), nothing ) end -function feynman_diagrams(proc::PROC, e::Int, u::Int, t::Int, m::Int) where {PROC<:AbstractProcessDefinition} - if (e + u + t) * 2 + m < 4 - throw(InvalidInputError("At least 4 particles must participate in a scattering process for a physical result!")) - end +function feynman_diagrams(proc::PROC) where {PROC<:GenericQEDProcess} + return feynman_diagrams(incoming_particles(proc), outgoing_particles(proc)) +end + +function feynman_diagrams(in_particles::Tuple, out_particles::Tuple) + _assert_particle_type_tuple(in_particles) + _assert_particle_type_tuple(out_particles) + + count(x -> x isa Electron, in_particles) + count(x -> x isa Positron, out_particles) == + count(x -> x isa Positron, in_particles) + count(x -> x isa Electron, out_particles) || + throw(InvalidInputError("the given particles do not make a physical process")) + + # get the fermion line counts and external photon count + e = count(x -> x isa Electron, in_particles) + count(x -> x isa Positron, out_particles) + m = count(x -> x isa Photon, in_particles) + count(x -> x isa Photon, out_particles) + # TODO: do this the same way as for e when muons and tauons are a part of QED.jl + u = 0 + t = 0 f_iter = _feynman_structures(e + u + t, m) e_perms = collect(permutations(Int[1:e;])) @@ -226,5 +234,5 @@ function feynman_diagrams(proc::PROC, e::Int, u::Int, t::Int, m::Int) where {PRO t_perms = collect(permutations(Int[e+u+1:e+u+t;])) first_photon_structure, _ = iterate(f_iter) - return FeynmanDiagramIterator(proc, Val(e), e_perms, 1, Val(u), u_perms, 1, Val(t), t_perms, 1, Val(m), f_iter, first_photon_structure) + return FeynmanDiagramIterator(Val(e), e_perms, 1, Val(u), u_perms, 1, Val(t), t_perms, 1, Val(m), f_iter, first_photon_structure) end diff --git a/src/generic_process_def.jl b/src/generic_process_def.jl new file mode 100644 index 0000000..311b1cc --- /dev/null +++ b/src/generic_process_def.jl @@ -0,0 +1,34 @@ +using QEDprocesses + +_assert_particle_type_tuple(::Tuple{}) = nothing +_assert_particle_type_tuple(t::Tuple{AbstractParticleType,Vararg}) = _assert_particle_type_tuple(t[2:end]) +_assert_particle_type_tuple(t::Any) = throw(InvalidInputError("invalid input, provide a tuple of AbstractParticleTypes to construct a GenericQEDProcess")) + +struct GenericQEDProcess{INT,OUTT} <: AbstractProcessDefinition where {INT<:Tuple,OUTT<:Tuple} + incoming_particles::INT + outgoing_particles::OUTT + + function GenericQEDProcess(in_particles::INT, out_particles::OUTT) where {INT<:Tuple,OUTT<:Tuple} + _assert_particle_type_tuple(in_particles) + _assert_particle_type_tuple(out_particles) + + #feynman_diagrams((), ()) + + return new{INT,OUTT}(in_particles, out_particles) + end +end + +QEDprocesses.incoming_particles(proc::GenericQEDProcess) = proc.incoming_particles +QEDprocesses.outgoing_particles(proc::GenericQEDProcess) = proc.outgoing_particles + +function isphysical(proc::GenericQEDProcess) + return (number_particles(proc, Incoming(), Electron()) + number_particles(proc, Outgoing(), Positron()) == + number_particles(proc, Incoming(), Positron()) + number_particles(proc, Outgoing(), Electron())) && + number_particles(proc, Incoming()) + number_particles(proc, Outgoing()) >= 2 +end + +function matrix_element(proc::GenericQEDProcess, psp::PhaseSpacePoint) + + + return nothing +end