From 70ce920ef69e242bbf55e780fb67d8fd5fe55641 Mon Sep 17 00:00:00 2001 From: Uwe Hernandez Acosta Date: Thu, 16 May 2024 15:53:49 +0200 Subject: [PATCH 01/22] updated tests for process interface --- test/interfaces/process_interface.jl | 85 ++++++++++-------------- test/test_implementation/test_process.jl | 71 +++++++++++--------- 2 files changed, 77 insertions(+), 79 deletions(-) diff --git a/test/interfaces/process_interface.jl b/test/interfaces/process_interface.jl index afdfa9e..b1a3d87 100644 --- a/test/interfaces/process_interface.jl +++ b/test/interfaces/process_interface.jl @@ -16,47 +16,54 @@ include("../test_implementation/TestImplementation.jl") OUTGOING_PARTICLES = rand(RNG, TestImplementation.PARTICLE_SET, N_OUTGOING) TESTPROC = TestImplementation.TestProcess(INCOMING_PARTICLES, OUTGOING_PARTICLES) - TESTPROC_FAIL = TestImplementation.TestProcess_FAIL( - INCOMING_PARTICLES, OUTGOING_PARTICLES - ) TESTMODEL = TestImplementation.TestModel() - TESTMODEL_FAIL = TestImplementation.TestModel_FAIL() TESTPSDEF = TestImplementation.TestPhasespaceDef() - TESTPSDEF_FAIL = TestImplementation.TestPhasespaceDef_FAIL() IN_PS = TestImplementation._rand_momenta(RNG, N_INCOMING) OUT_PS = TestImplementation._rand_momenta(RNG, N_OUTGOING) + PSP = generate_phase_space(TESTPROC,TESTMODEL,TESTPSDEF,IN_PS,OUT_PS) + @testset "failed interface" begin + TESTPROC_FAIL_ALL= TestImplementation.TestProcess_FAIL_ALL( + INCOMING_PARTICLES, OUTGOING_PARTICLES + ) + TESTPROC_FAIL_DIFFCS= TestImplementation.TestProcess_FAIL_DIFFCS( + INCOMING_PARTICLES, OUTGOING_PARTICLES + ) + TESTMODEL_FAIL = TestImplementation.TestModel_FAIL() + TESTPSDEF_FAIL = TestImplementation.TestPhasespaceDef_FAIL() + @testset "failed process interface" begin - @test_throws MethodError incoming_particles(TESTPROC_FAIL) - @test_throws MethodError outgoing_particles(TESTPROC_FAIL) + @test_throws MethodError incoming_particles(TESTPROC_FAIL_ALL) + @test_throws MethodError outgoing_particles(TESTPROC_FAIL_ALL) end + @testset "$PROC $MODEL" for (PROC, MODEL) in Iterators.product( - (TESTPROC, TESTPROC_FAIL), (TESTMODEL, TESTMODEL_FAIL) + (TESTPROC, TESTPROC_FAIL_DIFFCS), (TESTMODEL, TESTMODEL_FAIL) ) - in_ps = TestImplementation._rand_momenta(RNG, 2) - out_ps = TestImplementation._rand_momenta(RNG, 2) + if TestImplementation._any_fail(PROC, MODEL) - @test_throws MethodError QEDprocesses._incident_flux(PROC, MODEL, in_ps) - @test_throws MethodError QEDprocesses._averaging_norm(PROC, MODEL) - @test_throws MethodError QEDprocesses._matrix_element( - PROC, MODEL, in_ps, out_ps - ) + psp = generate_phase_space(PROC,MODEL,TESTPSDEF,IN_PS,OUT_PS) + @test_throws MethodError QEDprocesses._incident_flux(psp) + @test_throws MethodError QEDprocesses._averaging_norm(psp) + @test_throws MethodError QEDprocesses._matrix_element(psp) end for PS_DEF in (TESTPSDEF, TESTPSDEF_FAIL) if TestImplementation._any_fail(PROC, MODEL, PS_DEF) - @test_throws MethodError QEDprocesses._phase_space_factor( - PROC, MODEL, PS_DEF, in_ps, out_ps - ) + psp = generate_phase_space(PROC,MODEL,PS_DEF,IN_PS,OUT_PS) + @test_throws MethodError QEDprocesses._phase_space_factor(psp) end end end end @testset "broadcast" begin - test_func(proc) = proc + test_func(proc::AbstractProcessDefinition) = proc @test test_func.(TESTPROC) == TESTPROC + + test_func(model::AbstractModelDefinition) = model + @test test_func.(TESTMODEL) == TESTMODEL end @testset "incoming/outgoing particles" begin @@ -67,7 +74,7 @@ include("../test_implementation/TestImplementation.jl") end @testset "incident flux" begin - test_incident_flux = QEDprocesses._incident_flux(TESTPROC, TESTMODEL, IN_PS) + test_incident_flux = QEDprocesses._incident_flux(PSP) groundtruth = TestImplementation._groundtruth_incident_flux(IN_PS) @test isapprox(test_incident_flux, groundtruth, atol=ATOL, rtol=RTOL) end @@ -79,9 +86,7 @@ include("../test_implementation/TestImplementation.jl") end @testset "matrix element" begin - test_matrix_element = QEDprocesses._matrix_element( - TESTPROC, TESTMODEL, IN_PS, OUT_PS - ) + test_matrix_element = QEDprocesses._matrix_element(PSP) groundtruth = TestImplementation._groundtruth_matrix_element(IN_PS, OUT_PS) @test length(test_matrix_element) == length(groundtruth) for i in eachindex(test_matrix_element) @@ -90,41 +95,23 @@ include("../test_implementation/TestImplementation.jl") end @testset "is in phasespace" begin - @test QEDprocesses._is_in_phasespace(TESTPROC, TESTMODEL, TESTPSDEF, IN_PS, OUT_PS) - - IN_PS_unphysical = deepcopy(IN_PS) - IN_PS_unphysical[1] = SFourMomentum(zeros(4)) - - @test !QEDprocesses._is_in_phasespace( - TESTPROC, TESTMODEL, TESTPSDEF, IN_PS_unphysical, OUT_PS - ) - end - - @testset "is in phasespace" begin - @test QEDprocesses._is_in_phasespace(TESTPROC, TESTMODEL, TESTPSDEF, IN_PS, OUT_PS) + @test QEDprocesses._is_in_phasespace(PSP) IN_PS_unphysical = deepcopy(IN_PS) IN_PS_unphysical[1] = SFourMomentum(zeros(4)) - - @test !QEDprocesses._is_in_phasespace( - TESTPROC, TESTMODEL, TESTPSDEF, IN_PS_unphysical, OUT_PS - ) - OUT_PS_unphysical = deepcopy(OUT_PS) OUT_PS_unphysical[end] = ones(SFourMomentum) + PSP_unphysical_in_ps = generate_phase_space(TESTPROC,TESTMODEL,TESTPSDEF,IN_PS_unphysical,OUT_PS) + PSP_unphysical_out_ps = generate_phase_space(TESTPROC,TESTMODEL,TESTPSDEF,IN_PS,OUT_PS_unphysical) + PSP_unphysical= generate_phase_space(TESTPROC,TESTMODEL,TESTPSDEF,IN_PS_unphysical,OUT_PS_unphysical) - @test !QEDprocesses._is_in_phasespace( - TESTPROC, TESTMODEL, TESTPSDEF, IN_PS, OUT_PS_unphysical - ) - @test !QEDprocesses._is_in_phasespace( - TESTPROC, TESTMODEL, TESTPSDEF, IN_PS_unphysical, OUT_PS_unphysical - ) + @test !QEDprocesses._is_in_phasespace(PSP_unphysical_in_ps) + @test !QEDprocesses._is_in_phasespace(PSP_unphysical_out_ps) + @test !QEDprocesses._is_in_phasespace(PSP_unphysical) end @testset "phase space factor" begin - test_phase_space_factor = QEDprocesses._phase_space_factor( - TESTPROC, TESTMODEL, TESTPSDEF, IN_PS, OUT_PS - ) + test_phase_space_factor = QEDprocesses._phase_space_factor(PSP) groundtruth = TestImplementation._groundtruth_phase_space_factor(IN_PS, OUT_PS) @test isapprox(test_phase_space_factor, groundtruth, atol=ATOL, rtol=RTOL) end diff --git a/test/test_implementation/test_process.jl b/test/test_implementation/test_process.jl index fb7722d..f308681 100644 --- a/test/test_implementation/test_process.jl +++ b/test/test_implementation/test_process.jl @@ -24,33 +24,53 @@ end QEDprocesses.incoming_particles(proc::TestProcess) = proc.incoming_particles QEDprocesses.outgoing_particles(proc::TestProcess) = proc.outgoing_particles -struct TestProcess_FAIL{IP<:AbstractVector,OP<:AbstractVector} <: AbstractProcessDefinition +function QEDprocesses.in_phase_space_dimension(proc::TestProcess, ::TestModel) + return number_incoming_particles(proc) * 4 +end +function QEDprocesses.out_phase_space_dimension(proc::TestProcess, ::TestModel) + return number_outgoing_particles(proc) * 4 +end + +""" +Test process with no implemented interface. Should fail every usage except construction. +""" +struct TestProcess_FAIL_ALL{IP<:AbstractVector,OP<:AbstractVector} <: AbstractProcessDefinition incoming_particles::IP outgoing_particles::OP end -function TestProcess_FAIL(rng::AbstractRNG, N_in::Int, N_out::Int) +function TestProcess_FAIL_ALL(rng::AbstractRNG, N_in::Int, N_out::Int) in_particles = rand(rng, PARTICLE_SET, N_in) out_particles = rand(rng, PARTICLE_SET, N_out) - return TestProcess_FAIL(in_particles, out_particles) + return TestProcess_FAIL_ALL(in_particles, out_particles) end -function QEDprocesses.in_phase_space_dimension(proc::TestProcess, ::TestModel) - return number_incoming_particles(proc) * 4 +""" +Test process with no implemented interface except the incoming and outgoing particles. +Should fail every usage except construction of itself and the respective phase space point for given four-momenta. +""" +struct TestProcess_FAIL_DIFFCS{IP<:AbstractVector,OP<:AbstractVector} <: AbstractProcessDefinition + incoming_particles::IP + outgoing_particles::OP end -function QEDprocesses.out_phase_space_dimension(proc::TestProcess, ::TestModel) - return number_outgoing_particles(proc) * 4 + +function TestProcess_FAIL_DIFFCS(rng::AbstractRNG, N_in::Int, N_out::Int) + in_particles = rand(rng, PARTICLE_SET, N_in) + out_particles = rand(rng, PARTICLE_SET, N_out) + return TestProcess_FAIL_DIFFCS(in_particles, out_particles) end +QEDprocesses.incoming_particles(proc::TestProcess_FAIL_DIFFCS) = proc.incoming_particles +QEDprocesses.outgoing_particles(proc::TestProcess_FAIL_DIFFCS) = proc.outgoing_particles + # dummy phase space definition + failing phase space definition struct TestPhasespaceDef <: AbstractPhasespaceDefinition end struct TestPhasespaceDef_FAIL <: AbstractPhasespaceDefinition end # dummy implementation of the process interface -function QEDprocesses._incident_flux( - ::TestProcess, ::TestModel, in_ps::AbstractVector{T} -) where {T<:QEDbase.AbstractFourMomentum} +function QEDprocesses._incident_flux(psp::PhaseSpacePoint{<:TestProcess,TestModel}) + in_ps = momentum.(psp.in_particles) return _groundtruth_incident_flux(in_ps) end @@ -58,29 +78,21 @@ function QEDprocesses._averaging_norm(proc::TestProcess) return _groundtruth_averaging_norm(proc) end -function QEDprocesses._matrix_element( - ::TestProcess, ::TestModel, in_ps::AbstractVector{T}, out_ps::AbstractVector{T} -) where {T<:QEDbase.AbstractFourMomentum} +function QEDprocesses._matrix_element(psp::PhaseSpacePoint{<:TestProcess,TestModel}) + in_ps = momentum.(psp.in_particles) + out_ps = momentum.(psp.out_particles) return _groundtruth_matrix_element(in_ps, out_ps) end -function QEDprocesses._is_in_phasespace( - ::TestProcess, - ::TestModel, - ps_def::TestPhasespaceDef, - in_ps::AbstractVector{T}, - out_ps::AbstractVector{T}, -) where {T<:QEDbase.AbstractFourMomentum} +function QEDprocesses._is_in_phasespace(psp::PhaseSpacePoint{<:TestProcess,TestModel}) + in_ps = momentum.(psp.in_particles) + out_ps = momentum.(psp.out_particles) return _groundtruth_is_in_phasespace(in_ps, out_ps) end -function QEDprocesses._phase_space_factor( - ::TestProcess, - ::TestModel, - ps_def::TestPhasespaceDef, - in_ps::AbstractVector{T}, - out_ps::AbstractVector{T}, -) where {T<:QEDbase.AbstractFourMomentum} +function QEDprocesses._phase_space_factor(psp::PhaseSpacePoint{<:TestProcess,TestModel,TestPhasespaceDef}) + in_ps = momentum.(psp.in_particles) + out_ps = momentum.(psp.out_particles) return _groundtruth_phase_space_factor(in_ps, out_ps) end @@ -102,8 +114,7 @@ function QEDprocesses._generate_outgoing_momenta( return _groundtruth_generate_momenta(out_phase_space) end -function QEDprocesses._total_probability( - proc::TestProcess, model::TestModel, ps_def::TestPhasespaceDef, in_ps::AbstractVector{T} -) where {T<:QEDbase.AbstractFourMomentum} +function QEDprocesses._total_probability(psp::PhaseSpacePoint{<:TestProcess,TestModel}) + in_ps = momentum.(psp.in_particles) return _groundtruth_total_probability(in_ps) end From 6baf36c3968a9b48d74916a91725e54192d84b1e Mon Sep 17 00:00:00 2001 From: Uwe Hernandez Acosta Date: Thu, 16 May 2024 16:04:01 +0200 Subject: [PATCH 02/22] updated docstrings describing the process interface --- src/interfaces/process_interface.jl | 85 +++++++++-------------------- 1 file changed, 25 insertions(+), 60 deletions(-) diff --git a/src/interfaces/process_interface.jl b/src/interfaces/process_interface.jl index 8400ac4..70b5a9a 100644 --- a/src/interfaces/process_interface.jl +++ b/src/interfaces/process_interface.jl @@ -18,53 +18,26 @@ outgoing_particles(proc_def::AbstractProcessDefinition) which return a tuple of the incoming and outgoing particles, respectively. Furthermore, to calculate scattering probabilities and differential cross sections, the following -interface functions need to be implemented +interface functions need to be implemented for every combination of `CustomProcess<:AbstractProcessDefinition`, +`CustomModel<:AbstractModelDefinition`, and `CustomPhasespaceDefinition<:AbstractPhasespaceDefinition`. ```Julia - _incident_flux( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - in_phase_space::AbstractVector{T} - ) where {T<:QEDbase.AbstractFourMomentum} + _incident_flux(psp::PhaseSpacePoint{CustomProcess,CustomModel}) - _matrix_element( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - in_phase_space::AbstractVector{T} - out_phase_space::AbstractVector{T} - ) where {T<:QEDbase.AbstractFourMomentum} + _matrix_element(psp::PhaseSpacePoint{CustomProcess,CustomModel}) - _averaging_norm( - proc::AbstractProcessDefinition - ) + _averaging_norm(proc::CustomProcess) - _is_in_phasespace( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - ps_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVector{T} - out_phase_space::AbstractVector{T} - ) + _is_in_phasespace(psp::PhaseSpacePoint{CustomProcess,CustomModel}) - _phase_space_factor( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - ps_def::InPhasespaceDefinition, - in_phase_space::AbstractVector{T} - out_phase_space::AbstractVector{T} - ) where {T<:QEDbase.AbstractFourMomentum} + _phase_space_factor(psp::PhaseSpacePoint{CustomProcess,CustomModel, CustomPhasespaceDefinition}) ``` Optional is the implementation of ```Julia - _total_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - ps_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVector{T} - ) where {T<: QEDbase.AbstractFourMomentum} + _total_probability(psp::PhaseSpacePoint{CustomProcess,CustomModel, CustomPhasespaceDefinition}) ``` to enable the calculation of total probabilities and cross sections. @@ -94,11 +67,10 @@ This function needs to be given to implement the scattering process interface. function outgoing_particles end """ - _incident_flux( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - in_phase_space::AbstractVector{T} - ) where {T<:QEDbase.AbstractFourMomentum} + _incident_flux(PhaseSpacePoint{PROC,MODEL}) where { + PROC <: AbstractProcessDefinition, + MODEL <: AbstractModelDefinition + } Interface function, which returns the incident flux of the given scattering process for a given incoming phase-space. @@ -106,12 +78,10 @@ Interface function, which returns the incident flux of the given scattering proc function _incident_flux end """ - _matrix_element( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - in_phase_space::AbstractVector{T} - out_phase_space::AbstractVector{T} - ) where {T<:QEDbase.AbstractFourMomentum} + _matrix_elemen(PhaseSpacePoint{PROC,MODEL}) where { + PROC <: AbstractProcessDefinition, + MODEL <: AbstractModelDefinition, + } Interface function, which returns a tuple of scattering matrix elements for each spin and polarisation combination of `proc`. """ @@ -138,13 +108,10 @@ function _averaging_norm end """ - _is_in_phasespace( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - ps_def::AbstractPhasespaceDefinition, - in_ps::AbstractVector{T} - out_ps::AbstractVector{T} - ) where {T<:QEDbase.AbstractFourMomentum} + _is_in_phasespace(PhaseSpacePoint{PROC,MODEL}) where { + PROC <: AbstractProcessDefinition, + MODEL <: AbstractModelDefinition, + } Interface function, which returns `true`, if the combination of the given incoming and outgoing phase space is physical, i.e. all momenta are on-shell and some sort of energy-momentum conservation holds. @@ -152,13 +119,11 @@ is physical, i.e. all momenta are on-shell and some sort of energy-momentum cons function _is_in_phasespace end """ - _phase_space_factor( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - ps_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVector{T} - out_phase_space::AbstractVector{T} - ) where {T<:QEDbase.AbstractFourMomentum} + _phase_space_factor(PhaseSpacePoint{PROC,MODEL,PSDEF}) where { + PROC <: AbstractProcessDefinition, + MODEL <: AbstractModelDefinition + PSDEF <: AbstractPhasespaceDefinition, + } Interface function, which returns the pre-differential factor of the invariant phase space intergral measure. From ae0f18e6c8d16f9dbd5f8c7e1f3365b15a673056 Mon Sep 17 00:00:00 2001 From: Uwe Hernandez Acosta Date: Thu, 16 May 2024 16:08:56 +0200 Subject: [PATCH 03/22] formatting --- test/interfaces/process_interface.jl | 28 ++++++++++++++---------- test/test_implementation/test_process.jl | 12 ++++++---- 2 files changed, 24 insertions(+), 16 deletions(-) diff --git a/test/interfaces/process_interface.jl b/test/interfaces/process_interface.jl index b1a3d87..445ceb3 100644 --- a/test/interfaces/process_interface.jl +++ b/test/interfaces/process_interface.jl @@ -20,14 +20,13 @@ include("../test_implementation/TestImplementation.jl") TESTPSDEF = TestImplementation.TestPhasespaceDef() IN_PS = TestImplementation._rand_momenta(RNG, N_INCOMING) OUT_PS = TestImplementation._rand_momenta(RNG, N_OUTGOING) - PSP = generate_phase_space(TESTPROC,TESTMODEL,TESTPSDEF,IN_PS,OUT_PS) - + PSP = generate_phase_space(TESTPROC, TESTMODEL, TESTPSDEF, IN_PS, OUT_PS) @testset "failed interface" begin - TESTPROC_FAIL_ALL= TestImplementation.TestProcess_FAIL_ALL( + TESTPROC_FAIL_ALL = TestImplementation.TestProcess_FAIL_ALL( INCOMING_PARTICLES, OUTGOING_PARTICLES ) - TESTPROC_FAIL_DIFFCS= TestImplementation.TestProcess_FAIL_DIFFCS( + TESTPROC_FAIL_DIFFCS = TestImplementation.TestProcess_FAIL_DIFFCS( INCOMING_PARTICLES, OUTGOING_PARTICLES ) TESTMODEL_FAIL = TestImplementation.TestModel_FAIL() @@ -41,9 +40,8 @@ include("../test_implementation/TestImplementation.jl") @testset "$PROC $MODEL" for (PROC, MODEL) in Iterators.product( (TESTPROC, TESTPROC_FAIL_DIFFCS), (TESTMODEL, TESTMODEL_FAIL) ) - if TestImplementation._any_fail(PROC, MODEL) - psp = generate_phase_space(PROC,MODEL,TESTPSDEF,IN_PS,OUT_PS) + psp = generate_phase_space(PROC, MODEL, TESTPSDEF, IN_PS, OUT_PS) @test_throws MethodError QEDprocesses._incident_flux(psp) @test_throws MethodError QEDprocesses._averaging_norm(psp) @test_throws MethodError QEDprocesses._matrix_element(psp) @@ -51,7 +49,7 @@ include("../test_implementation/TestImplementation.jl") for PS_DEF in (TESTPSDEF, TESTPSDEF_FAIL) if TestImplementation._any_fail(PROC, MODEL, PS_DEF) - psp = generate_phase_space(PROC,MODEL,PS_DEF,IN_PS,OUT_PS) + psp = generate_phase_space(PROC, MODEL, PS_DEF, IN_PS, OUT_PS) @test_throws MethodError QEDprocesses._phase_space_factor(psp) end end @@ -63,7 +61,7 @@ include("../test_implementation/TestImplementation.jl") @test test_func.(TESTPROC) == TESTPROC test_func(model::AbstractModelDefinition) = model - @test test_func.(TESTMODEL) == TESTMODEL + @test test_func.(TESTMODEL) == TESTMODEL end @testset "incoming/outgoing particles" begin @@ -101,12 +99,18 @@ include("../test_implementation/TestImplementation.jl") IN_PS_unphysical[1] = SFourMomentum(zeros(4)) OUT_PS_unphysical = deepcopy(OUT_PS) OUT_PS_unphysical[end] = ones(SFourMomentum) - PSP_unphysical_in_ps = generate_phase_space(TESTPROC,TESTMODEL,TESTPSDEF,IN_PS_unphysical,OUT_PS) - PSP_unphysical_out_ps = generate_phase_space(TESTPROC,TESTMODEL,TESTPSDEF,IN_PS,OUT_PS_unphysical) - PSP_unphysical= generate_phase_space(TESTPROC,TESTMODEL,TESTPSDEF,IN_PS_unphysical,OUT_PS_unphysical) + PSP_unphysical_in_ps = generate_phase_space( + TESTPROC, TESTMODEL, TESTPSDEF, IN_PS_unphysical, OUT_PS + ) + PSP_unphysical_out_ps = generate_phase_space( + TESTPROC, TESTMODEL, TESTPSDEF, IN_PS, OUT_PS_unphysical + ) + PSP_unphysical = generate_phase_space( + TESTPROC, TESTMODEL, TESTPSDEF, IN_PS_unphysical, OUT_PS_unphysical + ) @test !QEDprocesses._is_in_phasespace(PSP_unphysical_in_ps) - @test !QEDprocesses._is_in_phasespace(PSP_unphysical_out_ps) + @test !QEDprocesses._is_in_phasespace(PSP_unphysical_out_ps) @test !QEDprocesses._is_in_phasespace(PSP_unphysical) end diff --git a/test/test_implementation/test_process.jl b/test/test_implementation/test_process.jl index f308681..30f9781 100644 --- a/test/test_implementation/test_process.jl +++ b/test/test_implementation/test_process.jl @@ -34,7 +34,8 @@ end """ Test process with no implemented interface. Should fail every usage except construction. """ -struct TestProcess_FAIL_ALL{IP<:AbstractVector,OP<:AbstractVector} <: AbstractProcessDefinition +struct TestProcess_FAIL_ALL{IP<:AbstractVector,OP<:AbstractVector} <: + AbstractProcessDefinition incoming_particles::IP outgoing_particles::OP end @@ -49,7 +50,8 @@ end Test process with no implemented interface except the incoming and outgoing particles. Should fail every usage except construction of itself and the respective phase space point for given four-momenta. """ -struct TestProcess_FAIL_DIFFCS{IP<:AbstractVector,OP<:AbstractVector} <: AbstractProcessDefinition +struct TestProcess_FAIL_DIFFCS{IP<:AbstractVector,OP<:AbstractVector} <: + AbstractProcessDefinition incoming_particles::IP outgoing_particles::OP end @@ -69,7 +71,7 @@ struct TestPhasespaceDef_FAIL <: AbstractPhasespaceDefinition end # dummy implementation of the process interface -function QEDprocesses._incident_flux(psp::PhaseSpacePoint{<:TestProcess,TestModel}) +function QEDprocesses._incident_flux(psp::PhaseSpacePoint{<:TestProcess,TestModel}) in_ps = momentum.(psp.in_particles) return _groundtruth_incident_flux(in_ps) end @@ -90,7 +92,9 @@ function QEDprocesses._is_in_phasespace(psp::PhaseSpacePoint{<:TestProcess,TestM return _groundtruth_is_in_phasespace(in_ps, out_ps) end -function QEDprocesses._phase_space_factor(psp::PhaseSpacePoint{<:TestProcess,TestModel,TestPhasespaceDef}) +function QEDprocesses._phase_space_factor( + psp::PhaseSpacePoint{<:TestProcess,TestModel,TestPhasespaceDef} +) in_ps = momentum.(psp.in_particles) out_ps = momentum.(psp.out_particles) return _groundtruth_phase_space_factor(in_ps, out_ps) From f9705ae5bb29acc9dfb19ec72ad069d9de7d131c Mon Sep 17 00:00:00 2001 From: Uwe Hernandez Acosta Date: Thu, 16 May 2024 22:18:27 +0200 Subject: [PATCH 04/22] opt out the evaluation of cross sections and probabilities on sets of inputs --- src/cross_sections.jl | 70 ---------------------------------------- src/probabilities.jl | 72 ------------------------------------------ test/cross_sections.jl | 63 ++++++------------------------------ test/runtests.jl | 2 +- 4 files changed, 10 insertions(+), 197 deletions(-) diff --git a/src/cross_sections.jl b/src/cross_sections.jl index bf3279a..d149d52 100644 --- a/src/cross_sections.jl +++ b/src/cross_sections.jl @@ -44,42 +44,6 @@ function _unsafe_differential_cross_section( ) end -# differential cross sections without energy momentum conservation check -# single in phase space point/ several out phase space points -function _unsafe_differential_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVector{T}, - out_phase_space::AbstractMatrix{T}, -) where {T<:AbstractPhasespaceElement} - res = Vector{eltype(T)}(undef, size(out_phase_space, 2)) - for i in 1:size(out_phase_space, 2) - res[i] = _unsafe_differential_cross_section( - proc, model, phase_space_def, in_phase_space, view(out_phase_space, :, i) - ) - end - return res -end - -# differential cross sections without energy momentum conservation check -# several in phase space points/ one or several out phase space points -function _unsafe_differential_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractMatrix{T}, - out_phase_space::AbstractVecOrMat{T}, -) where {T<:AbstractPhasespaceElement} - res = Matrix{eltype(T)}(undef, size(in_phase_space, 2), size(out_phase_space, 2)) - for i in 1:size(in_phase_space, 2) - res[i, :] .= _unsafe_differential_cross_section( - proc, model, phase_space_def, view(in_phase_space, :, i), out_phase_space - ) - end - return res -end - """ function unsafe_differential_cross_section( @@ -183,40 +147,6 @@ function _differential_cross_section( ) end -function _differential_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVector{T}, - out_phase_space::AbstractMatrix{T}, -) where {T<:AbstractPhasespaceElement} - res = Vector{eltype(T)}(undef, size(out_phase_space, 2)) - for i in 1:size(out_phase_space, 2) - res[i] = _differential_cross_section( - proc, model, phase_space_def, in_phase_space, view(out_phase_space, :, i) - ) - end - return res -end - -# differential cross sections with energy momentum conservation check -# several in phase space points/ one or several out phase space points -function _differential_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractMatrix{T}, - out_phase_space::AbstractVecOrMat{T}, -) where {T<:AbstractPhasespaceElement} - res = Matrix{eltype(T)}(undef, size(in_phase_space, 2), size(out_phase_space, 2)) - for i in 1:size(in_phase_space, 2) - res[i, :] .= _differential_cross_section( - proc, model, phase_space_def, view(in_phase_space, :, i), out_phase_space - ) - end - return res -end - """ differential_cross_section( proc::AbstractProcessDefinition, diff --git a/src/probabilities.jl b/src/probabilities.jl index b5f24e0..cee0181 100644 --- a/src/probabilities.jl +++ b/src/probabilities.jl @@ -46,42 +46,6 @@ function _unsafe_differential_probability( ) end -# differential probability without energy momentum conservation check -# single in phase space points/ several out phase space point -function _unsafe_differential_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVector{T}, - out_phase_space::AbstractMatrix{T}, -) where {T<:AbstractPhasespaceElement} - res = Vector{eltype(T)}(undef, size(out_phase_space, 2)) - for i in 1:size(out_phase_space, 2) - res[i] = _unsafe_differential_probability( - proc, model, phase_space_def, in_phase_space, view(out_phase_space, :, i) - ) - end - return res -end - -# differential probability without energy momentum conservation check -# several in phase space points/ one or several out phase space point -function _unsafe_differential_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractMatrix{T}, - out_phase_space::AbstractVecOrMat{T}, -) where {T<:AbstractPhasespaceElement} - res = Matrix{eltype(T)}(undef, size(in_phase_space, 2), size(out_phase_space, 2)) - for i in 1:size(in_phase_space, 2) - res[i, :] .= _unsafe_differential_probability( - proc, model, phase_space_def, view(in_phase_space, :, i), out_phase_space - ) - end - return res -end - """ unsafe_differential_probability(phase_space_point::PhaseSpacePoint) @@ -184,42 +148,6 @@ function _differential_probability( return _differential_probability(proc, model, phase_space_def, in_momenta, out_momenta) end -# differential probability with energy momentum conservation check -# one in phase space points/ several out phase space point -function _differential_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVector{T}, - out_phase_space::AbstractMatrix{T}, -) where {T<:AbstractPhasespaceElement} - res = Vector{eltype(T)}(undef, size(out_phase_space, 2)) - for i in 1:size(out_phase_space, 2) - res[i] = _differential_probability( - proc, model, phase_space_def, in_phase_space, view(out_phase_space, :, i) - ) - end - return res -end - -# differential probability with energy momentum conservation check -# several in phase space points/ one or several out phase space point -function _differential_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractMatrix{T}, - out_phase_space::AbstractVecOrMat{T}, -) where {T<:AbstractPhasespaceElement} - res = Matrix{eltype(T)}(undef, size(in_phase_space, 2), size(out_phase_space, 2)) - for i in 1:size(in_phase_space, 2) - res[i, :] .= _differential_probability( - proc, model, phase_space_def, view(in_phase_space, :, i), out_phase_space - ) - end - return res -end - """ differential_probability(phase_space_point::PhaseSpacePoint) diff --git a/test/cross_sections.jl b/test/cross_sections.jl index 42df677..cbc1423 100644 --- a/test/cross_sections.jl +++ b/test/cross_sections.jl @@ -30,48 +30,12 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() p_out_unphys = TestImplementation._rand_out_momenta_failing(RNG, N_OUTGOING) p_out_unphys_invalid = TestImplementation._rand_out_momenta_failing(RNG, N_OUTGOING + 1) - # sets of ps points - p_in_set_phys = TestImplementation._rand_momenta(RNG, N_INCOMING, 2) - p_in_set_unphys_mix = TestImplementation._rand_in_momenta_failing_mix( - RNG, N_INCOMING, 2 - ) - p_in_set_unphys_all = TestImplementation._rand_in_momenta_failing_all( - RNG, N_INCOMING, 2 - ) - p_in_set_phys_invalid = TestImplementation._rand_momenta(RNG, N_INCOMING + 1, 2) - p_in_set_unphys_mix_invalid = TestImplementation._rand_in_momenta_failing_mix( - RNG, N_INCOMING + 1, 2 - ) - p_in_set_unphys_all_invalid = TestImplementation._rand_in_momenta_failing_all( - RNG, N_INCOMING + 1, 2 - ) - - p_out_set_phys = TestImplementation._rand_momenta(RNG, N_OUTGOING, 2) - p_out_set_unphys_mix = TestImplementation._rand_out_momenta_failing_mix( - RNG, N_OUTGOING, 2 - ) - p_out_set_unphys_all = TestImplementation._rand_out_momenta_failing_all( - RNG, N_OUTGOING, 2 - ) - p_out_set_phys_invalid = TestImplementation._rand_momenta(RNG, N_OUTGOING + 1, 2) - p_out_set_unphys_mix_invalid = TestImplementation._rand_out_momenta_failing_mix( - RNG, N_OUTGOING + 1, 2 - ) - p_out_set_unphys_all_invalid = TestImplementation._rand_out_momenta_failing_all( - RNG, N_OUTGOING + 1, 2 - ) p_in_all = ( p_in_phys, p_in_unphys, p_in_phys_invalid, p_in_unphys_invalid, - p_in_set_phys, - p_in_set_unphys_mix, - p_in_set_unphys_all, - p_in_set_phys_invalid, - p_in_set_unphys_mix_invalid, - p_in_set_unphys_all_invalid, ) p_out_all = ( @@ -79,29 +43,20 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() p_out_phys_invalid, p_out_unphys, p_out_unphys_invalid, - p_out_set_phys, - p_out_set_unphys_mix, - p_out_set_unphys_all, - p_out_set_phys_invalid, - p_out_set_unphys_mix_invalid, - p_out_set_unphys_all_invalid, ) + # all combinations p_combs = Iterators.product(p_in_all, p_out_all) - p_in_all_valid = ( - p_in_phys, p_in_unphys, p_in_set_phys, p_in_set_unphys_mix, p_in_set_unphys_all - ) + p_in_all_valid = (p_in_phys, p_in_unphys) - p_out_all_valid = ( - p_out_phys, p_out_unphys, p_out_set_phys, p_out_set_unphys_mix, p_out_set_unphys_all - ) + p_out_all_valid = (p_out_phys, p_out_unphys) # all valid combinations p_combs_valid = Iterators.product(p_in_all_valid, p_out_all_valid) - p_in_all_phys = (p_in_phys, p_in_set_phys) - p_out_all_phys = (p_out_phys, p_out_set_phys) + p_in_all_phys = (p_in_phys, ) + p_out_all_phys = (p_out_phys, ) p_combs_phys = Iterators.product(p_in_all_phys, p_out_all_phys) @@ -200,7 +155,7 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() @testset "total cross section" begin @testset "compute" begin - for P_IN in (p_in_phys, p_in_set_phys) + for P_IN in (p_in_phys,) groundtruth = TestImplementation._groundtruth_total_cross_section( P_IN ) @@ -218,7 +173,7 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() end end @testset "invalid input" begin - for P_IN in (p_in_phys_invalid, p_in_set_phys_invalid) + for P_IN in (p_in_phys_invalid, ) @test_throws DimensionMismatch total_cross_section( TESTPROC, TESTMODEL, TESTPSDEF, P_IN ) @@ -333,7 +288,7 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() @testset "total probability" begin @testset "compute" begin - for P_IN in (p_in_phys, p_in_set_phys) + for P_IN in (p_in_phys,) groundtruth = TestImplementation._groundtruth_total_probability( P_IN ) @@ -351,7 +306,7 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() end end @testset "invalid input" begin - for P_IN in (p_in_phys_invalid, p_in_set_phys_invalid) + for P_IN in (p_in_phys_invalid) @test_throws DimensionMismatch total_probability( TESTPROC, TESTMODEL, TESTPSDEF, P_IN ) diff --git a/test/runtests.jl b/test/runtests.jl index 9ecc6ee..76c49ae 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -19,7 +19,7 @@ begin include("propagators.jl") end @time @safetestset "cross section & probability" begin - include("cross_sections.jl") + #include("cross_sections.jl") end @time @safetestset "phase spaces" begin From 24e49b23498b2a74e67d77cc4959176dfde7ec8a Mon Sep 17 00:00:00 2001 From: Uwe Hernandez Acosta Date: Thu, 16 May 2024 22:18:44 +0200 Subject: [PATCH 05/22] formatting --- test/cross_sections.jl | 21 +++++---------------- 1 file changed, 5 insertions(+), 16 deletions(-) diff --git a/test/cross_sections.jl b/test/cross_sections.jl index cbc1423..79d832d 100644 --- a/test/cross_sections.jl +++ b/test/cross_sections.jl @@ -30,20 +30,9 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() p_out_unphys = TestImplementation._rand_out_momenta_failing(RNG, N_OUTGOING) p_out_unphys_invalid = TestImplementation._rand_out_momenta_failing(RNG, N_OUTGOING + 1) + p_in_all = (p_in_phys, p_in_unphys, p_in_phys_invalid, p_in_unphys_invalid) - p_in_all = ( - p_in_phys, - p_in_unphys, - p_in_phys_invalid, - p_in_unphys_invalid, - ) - - p_out_all = ( - p_out_phys, - p_out_phys_invalid, - p_out_unphys, - p_out_unphys_invalid, - ) + p_out_all = (p_out_phys, p_out_phys_invalid, p_out_unphys, p_out_unphys_invalid) # all combinations p_combs = Iterators.product(p_in_all, p_out_all) @@ -55,8 +44,8 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() # all valid combinations p_combs_valid = Iterators.product(p_in_all_valid, p_out_all_valid) - p_in_all_phys = (p_in_phys, ) - p_out_all_phys = (p_out_phys, ) + p_in_all_phys = (p_in_phys,) + p_out_all_phys = (p_out_phys,) p_combs_phys = Iterators.product(p_in_all_phys, p_out_all_phys) @@ -173,7 +162,7 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() end end @testset "invalid input" begin - for P_IN in (p_in_phys_invalid, ) + for P_IN in (p_in_phys_invalid,) @test_throws DimensionMismatch total_cross_section( TESTPROC, TESTMODEL, TESTPSDEF, P_IN ) From 69593f204565202842e2eb8308291135ab78dcb2 Mon Sep 17 00:00:00 2001 From: Uwe Hernandez Acosta Date: Thu, 16 May 2024 23:34:36 +0200 Subject: [PATCH 06/22] updated diffCS and diffProb to PSPs, removed old interface, updated tests, rearranged files --- src/QEDprocesses.jl | 8 +- src/cross_section/diff_cross_section.jl | 43 ++++ src/cross_section/diff_probability.jl | 44 ++++ src/cross_section/internal.jl | 131 ++++++++++ src/cross_section/total_cross_section.jl | 87 +++++++ src/cross_section/total_probability.jl | 74 ++++++ src/cross_sections.jl | 303 ----------------------- src/interfaces/process_interface.jl | 24 +- src/probabilities.jl | 290 ---------------------- test/cross_sections.jl | 241 ++++-------------- test/interfaces/process_interface.jl | 2 +- test/test_implementation/test_process.jl | 10 +- 12 files changed, 453 insertions(+), 804 deletions(-) create mode 100644 src/cross_section/diff_cross_section.jl create mode 100644 src/cross_section/diff_probability.jl create mode 100644 src/cross_section/internal.jl create mode 100644 src/cross_section/total_cross_section.jl create mode 100644 src/cross_section/total_probability.jl delete mode 100644 src/cross_sections.jl delete mode 100644 src/probabilities.jl diff --git a/src/QEDprocesses.jl b/src/QEDprocesses.jl index d571b0f..755a8e3 100644 --- a/src/QEDprocesses.jl +++ b/src/QEDprocesses.jl @@ -58,8 +58,12 @@ include("interfaces/setup_interface.jl") include("phase_spaces.jl") include("momentum_generation.jl") include("propagators.jl") -include("probabilities.jl") -include("cross_sections.jl") + +include("cross_section/diff_probability.jl") +include("cross_section/diff_cross_section.jl") +include("cross_section/total_probability.jl") +include("cross_section/total_cross_section.jl") +include("cross_section/internal.jl") include("models/models.jl") include("processes/one_photon_compton/one_photon_compton.jl") diff --git a/src/cross_section/diff_cross_section.jl b/src/cross_section/diff_cross_section.jl new file mode 100644 index 0000000..74bcb8d --- /dev/null +++ b/src/cross_section/diff_cross_section.jl @@ -0,0 +1,43 @@ +######################## +# differential and total cross sections. +# +# This file contains default implementations for differential +# cross sections based on the scattering process interface +######################## + + +function _incident_flux(psp::PhaseSpacePoint) + return _incident_flux( + psp.proc, + psp.model, + momentum.(psp.in_particles) + ) +end + +""" + unsafe_differential_cross_section(phase_space_point::PhaseSpacePoint) + +Return the differential cross section evaluated on a phase space point without checking if the given phase space is physical. +""" +function unsafe_differential_cross_section(phase_space_point::PhaseSpacePoint) + I = 1 / (4 * _incident_flux(phase_space_point)) + + return I * unsafe_differential_probability(phase_space_point) +end + + +""" + differential_cross_section(phase_space_point::PhaseSpacePoint) + +If the given phase spaces are physical, return differential cross section evaluated on a phase space point. Zero otherwise. +""" +function differential_cross_section(phase_space_point::PhaseSpacePoint) + if !_is_in_phasespace(phase_space_point) + return zero(eltype(momentum(phase_space_point,Incoming(),1))) + end + + return unsafe_differential_cross_section(phase_space_point) +end + + + diff --git a/src/cross_section/diff_probability.jl b/src/cross_section/diff_probability.jl new file mode 100644 index 0000000..423f865 --- /dev/null +++ b/src/cross_section/diff_probability.jl @@ -0,0 +1,44 @@ +############ +# scattering probabilities +# +# This file contains implementations of the scattering probability based on the +# process interface with and without input validation and/or phase space +# constraint. +############ + +# convenience function +# can be overloaded if an analytical version is known +function _matrix_element_square(psp::PhaseSpacePoint) + mat_el = _matrix_element(psp) + return abs2.(mat_el) +end + +""" + unsafe_differential_probability(phase_space_point::PhaseSpacePoint) + +Return differential probability evaluated on a phase space point without checking if the given phase space(s) are physical. +""" +function unsafe_differential_probability(psp::PhaseSpacePoint) + matrix_elements_sq = _matrix_element_square(psp) + + normalization = _averaging_norm(psp.proc) + + ps_fac = _phase_space_factor(psp) + + return normalization * sum(matrix_elements_sq) * ps_fac +end + +""" + differential_probability(phase_space_point::PhaseSpacePoint) + +If the given phase spaces are physical, return differential probability evaluated on a phase space point. Zero otherwise. +""" +function differential_probability(phase_space_point::PhaseSpacePoint) + if !_is_in_phasespace(phase_space_point) + return zero(eltype(momentum(phase_space_point,Incoming(),1))) + end + + return unsafe_differential_probability(phase_space_point) +end + + diff --git a/src/cross_section/internal.jl b/src/cross_section/internal.jl new file mode 100644 index 0000000..d0e8f96 --- /dev/null +++ b/src/cross_section/internal.jl @@ -0,0 +1,131 @@ + +# convenience function for internal use only +# differential probability without energy momentum conservation check +# single in phase space points/ single out phase space point +# based on four-momenta +function _unsafe_differential_probability( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + phase_space_def::AbstractPhasespaceDefinition, + in_phase_space::AbstractVector{T}, + out_phase_space::AbstractVector{T}, +) where {T<:QEDbase.AbstractFourMomentum} + psp = generate_phase_space(proc,model,phase_space_def,in_phase_space,out_phase_space) + return unsafe_differential_probability(psp) +end + +# convenience function for internal use only +# differential probability without energy momentum conservation check +# based on coordinates +function _unsafe_differential_probability( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + phase_space_def::AbstractPhasespaceDefinition, + in_phase_space::AbstractVector{T}, + out_phase_space::AbstractVector{T}, +) where {T<:Real} + in_momenta, out_momenta = _generate_momenta( + proc, model, phase_space_def, in_phase_space, out_phase_space + ) + return _unsafe_differential_probability( + proc, model, phase_space_def, in_momenta, out_momenta + ) +end + +# convenience function for internal use only +# differential probability with energy momentum conservation check +# one in phase space point/ one out phase space point +# based on four-momenta +function _differential_probability( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + phase_space_def::AbstractPhasespaceDefinition, + in_phase_space::AbstractVector{T}, + out_phase_space::AbstractVector{T}, +) where {T<:QEDbase.AbstractFourMomentum} + psp = generate_phase_space(proc,model,phase_space_def,in_phase_space,out_phase_space) + return differential_probability(psp) +end + +# convenience function for internal use only +# differential probability with energy momentum conservation check +# one in phase space point/ one out phase space point +# based on coordinates +function _differential_probability( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + phase_space_def::AbstractPhasespaceDefinition, + in_phase_space::AbstractVector{T}, + out_phase_space::AbstractVector{T}, +) where {T<:Real} + in_momenta, out_momenta = _generate_momenta( + proc, model, phase_space_def, in_phase_space, out_phase_space + ) + return _differential_probability(proc, model, phase_space_def, in_momenta, out_momenta) +end + +# convenience function for internal use only +# differential cross sections without energy momentum conservation check +# single in phase space point/ single out phase space point +# based on four-momenta +function _unsafe_differential_cross_section( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + phase_space_def::AbstractPhasespaceDefinition, + in_phase_space::AbstractVector{T}, + out_phase_space::AbstractVector{T}, +) where {T<:QEDbase.AbstractFourMomentum} + psp = generate_phase_space(proc,model,phase_space_def,in_phase_space,out_phase_space) + return _unsafe_differential_cross_section(phase_space_point) +end + +# convenience function for internal use only +# differential cross sections without energy momentum conservation check +# single in phase space point/ single out phase space point +# based on coordinates +function _unsafe_differential_cross_section( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + phase_space_def::AbstractPhasespaceDefinition, + in_phase_space::AbstractVector{T}, + out_phase_space::AbstractVector{T}, +) where {T<:Real} + in_momenta, out_momenta = _generate_momenta( + proc, model, phase_space_def, in_phase_space, out_phase_space + ) + return _unsafe_differential_cross_section( + proc, model, phase_space_def, in_momenta, out_momenta + ) +end + +# convenience function for internal use only +# differential cross sections with energy momentum conservation check +# single in phase space point/ single out phase space point +function _differential_cross_section( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + phase_space_def::AbstractPhasespaceDefinition, + in_phase_space::AbstractVector{T}, + out_phase_space::AbstractVector{T}, +) where {T<:QEDbase.AbstractFourMomentum} + psp = generate_phase_space(proc,model,phase_space_def,in_phase_space,out_phase_space) + return differential_cross_section(psp) +end + +# convenience function for internal use only +# differential cross sections with energy momentum conservation check +# single in phase space point/ several out phase space points +function _differential_cross_section( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + phase_space_def::AbstractPhasespaceDefinition, + in_phase_space::AbstractVector{T}, + out_phase_space::AbstractVector{T}, +)::Float64 where {T<:Real} + in_momenta, out_momenta = _generate_momenta( + proc, model, phase_space_def, in_phase_space, out_phase_space + ) + return _differential_cross_section( + proc, model, phase_space_def, in_momenta, out_momenta + ) +end diff --git a/src/cross_section/total_cross_section.jl b/src/cross_section/total_cross_section.jl new file mode 100644 index 0000000..a7574dd --- /dev/null +++ b/src/cross_section/total_cross_section.jl @@ -0,0 +1,87 @@ + +############ +# Total cross sections +############ + +# total cross section on single phase space point +# based on four-momenta +function _total_cross_section( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + phase_space_def::AbstractPhasespaceDefinition, + in_phase_space::AbstractVector{T}, +) where {T<:QEDbase.AbstractFourMomentum} + I = 1 / (4 * _incident_flux(proc, model, in_phase_space)) + + return I * _total_probability(proc, model, phase_space_def, in_phase_space) +end + +# total cross section on single phase space point +# based on coordinates +function _total_cross_section( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + phase_space_def::AbstractPhasespaceDefinition, + in_phase_space::AbstractVector{T}, +) where {T<:Real} + in_momenta = _generate_incoming_momenta(proc, model, phase_space_def, in_phase_space) + return _total_cross_section(proc, model, phase_space_def, in_momenta) +end + +# total cross section on several phase space points +function _total_cross_section( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + in_phase_space_def::AbstractPhasespaceDefinition, + in_phase_space::AbstractMatrix{T}, +) where {T<:AbstractPhasespaceElement} + res = Vector{eltype(T)}(undef, size(in_phase_space, 2)) + for i in 1:size(in_phase_space, 2) + res[i] = _total_cross_section( + proc, model, in_phase_space_def, view(in_phase_space, :, i) + ) + end + return res +end + +""" + total_cross_section( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + in_phase_space_def::AbstractPhasespaceDefinition, + in_phase_space::AbstractVecOrMat{T}, + ) where {T<:QEDbase.AbstractFourMomentum} + +Return the total cross section for a given combination of scattering process and compute model, evaluated at the particle momenta. +""" +function total_cross_section( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + in_phase_space_def::AbstractPhasespaceDefinition, + in_phase_space::AbstractVecOrMat{T}, +) where {T<:QEDbase.AbstractFourMomentum} + _check_in_phase_space_dimension(proc, model, in_phase_space) + + return _total_cross_section(proc, model, in_phase_space_def, in_phase_space) +end + +""" + total_cross_section( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + in_phase_space_def::AbstractPhasespaceDefinition, + in_phase_space::AbstractVecOrMat{T}, + ) where {T<:Real} + +Return the total cross section for a given combination of scattering process and compute model, evaluated at the coordinates. +""" +function total_cross_section( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + in_phase_space_def::AbstractPhasespaceDefinition, + in_phase_space::AbstractVecOrMat{T}, +) where {T<:Real} + _check_in_phase_space_dimension(proc, model, in_phase_space) + + return _total_cross_section(proc, model, in_phase_space_def, in_phase_space) +end diff --git a/src/cross_section/total_probability.jl b/src/cross_section/total_probability.jl new file mode 100644 index 0000000..3296ba0 --- /dev/null +++ b/src/cross_section/total_probability.jl @@ -0,0 +1,74 @@ + +########### +# Total probability +########### + +# total probability on a phase space point +# based on coordinates +function _total_probability( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + phase_space_def::AbstractPhasespaceDefinition, + in_phase_space::AbstractVector{T}, +) where {T<:Real} + in_momenta = _generate_incoming_momenta(proc, model, phase_space_def, in_phase_space) + return _total_probability(proc, model, phase_space_def, in_momenta) +end + +# total probability on several phase space points +function _total_probability( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + phase_space_def::AbstractPhasespaceDefinition, + in_phase_space::AbstractMatrix{T}, +) where {T<:AbstractPhasespaceElement} + res = Vector{eltype(T)}(undef, size(in_phase_space, 2)) + for i in 1:size(in_phase_space, 2) + res[i] = _total_probability( + proc, model, phase_space_def, view(in_phase_space, :, i) + ) + end + return res +end + +""" + total_probability( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + phase_space_def::AbstractPhasespaceDefinition, + in_phase_space::AbstractMatrix{T}, + ) where {T<:QEDbase.AbstractFourMomentum} + +Return the total probability of a given model and process combination, evaluated at the particle momenta. +""" +function total_probability( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + phase_space_def::AbstractPhasespaceDefinition, + in_phase_space::AbstractVecOrMat{T}, +) where {T<:QEDbase.AbstractFourMomentum} + _check_in_phase_space_dimension(proc, model, in_phase_space) + + return _total_probability(proc, model, phase_space_def, in_phase_space) +end + +""" + total_probability( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + phase_space_def::AbstractPhasespaceDefinition, + in_phase_space::AbstractMatrix{T}, + ) where {T<:Real} + +Return the total probability of a given model and process combination, evaluated at the coordinates. +""" +function total_probability( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + phase_space_def::AbstractPhasespaceDefinition, + in_phase_space::AbstractVecOrMat{T}, +) where {T<:Real} + _check_in_phase_space_dimension(proc, model, in_phase_space) + + return _total_probability(proc, model, phase_space_def, in_phase_space) +end diff --git a/src/cross_sections.jl b/src/cross_sections.jl deleted file mode 100644 index d149d52..0000000 --- a/src/cross_sections.jl +++ /dev/null @@ -1,303 +0,0 @@ -######################## -# differential and total cross sections. -# -# This file contains default implementations for differential and total cross -# sections based on the scattering process interface -######################## - -############ -# differential cross sections -############ - -# differential cross sections without energy momentum conservation check -# single in phase space point/ single out phase space point -# based on four-momenta -function _unsafe_differential_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVector{T}, - out_phase_space::AbstractVector{T}, -) where {T<:QEDbase.AbstractFourMomentum} - I = 1 / (4 * _incident_flux(proc, model, in_phase_space)) - - return I * _unsafe_differential_probability( - proc, model, phase_space_def, in_phase_space, out_phase_space - ) -end - -# differential cross sections without energy momentum conservation check -# single in phase space point/ single out phase space point -# based on coordinates -function _unsafe_differential_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVector{T}, - out_phase_space::AbstractVector{T}, -) where {T<:Real} - in_momenta, out_momenta = _generate_momenta( - proc, model, phase_space_def, in_phase_space, out_phase_space - ) - return _unsafe_differential_cross_section( - proc, model, phase_space_def, in_momenta, out_momenta - ) -end - -""" - - function unsafe_differential_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, - out_phase_space::AbstractVecOrMat{T}, - ) where {T<:QEDbase.AbstractFourMomentum} - -Return the differential cross section evaluated at the four-momenta without checking if the given phase space(s) are physical. -""" -function unsafe_differential_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, - out_phase_space::AbstractVecOrMat{T}, -) where {T<:QEDbase.AbstractFourMomentum} - _check_in_phase_space_dimension(proc, model, in_phase_space) - _check_out_phase_space_dimension(proc, model, out_phase_space) - - return _unsafe_differential_cross_section( - proc, model, phase_space_def, in_phase_space, out_phase_space - ) -end - -""" - unsafe_differential_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, - out_phase_space::AbstractVecOrMat{T}, -) where {T<:Real} - -Return the differential cross section evaluated at the coordinates without checking if the given phase space(s) are physical. -""" -function unsafe_differential_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, - out_phase_space::AbstractVecOrMat{T}, -) where {T<:Real} - _check_in_phase_space_dimension(proc, model, in_phase_space) - _check_out_phase_space_dimension(proc, model, out_phase_space) - - return _unsafe_differential_cross_section( - proc, model, phase_space_def, in_phase_space, out_phase_space - ) -end - -""" - unsafe_differential_cross_section(phase_space_point::PhaseSpacePoint) - -Return the differential cross section evaluated on a phase space point without checking if the given phase space is physical. -""" -function unsafe_differential_cross_section(phase_space_point::PhaseSpacePoint) - return _unsafe_differential_cross_section( - phase_space_point.proc, - phase_space_point.model, - phase_space_point.ps_def, - momentum.(phase_space_point.in_particles), - momentum.(phase_space_point.out_particles), - ) -end - -# differential cross sections with energy momentum conservation check -# single in phase space point/ single out phase space point -function _differential_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVector{T}, - out_phase_space::AbstractVector{T}, -) where {T<:QEDbase.AbstractFourMomentum} - if !_is_in_phasespace(proc, model, phase_space_def, in_phase_space, out_phase_space) - return zero(eltype(T)) - end - - return _unsafe_differential_cross_section( - proc, model, phase_space_def, in_phase_space, out_phase_space - ) -end - -# differential cross sections with energy momentum conservation check -# single in phase space point/ several out phase space points -function _differential_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVector{T}, - out_phase_space::AbstractVector{T}, -)::Float64 where {T<:Real} - in_momenta, out_momenta = _generate_momenta( - proc, model, phase_space_def, in_phase_space, out_phase_space - ) - return _differential_cross_section( - proc, model, phase_space_def, in_momenta, out_momenta - ) -end - -""" - differential_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, - out_phase_space::AbstractVecOrMat{T}, - ) where {T<:QEDbase.AbstractFourMomentum} - -If the given phase spaces are physical, return differential cross section evaluated at the four-momenta. Zero otherwise. - -""" -function differential_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, - out_phase_space::AbstractVecOrMat{T}, -) where {T<:QEDbase.AbstractFourMomentum} - _check_in_phase_space_dimension(proc, model, in_phase_space) - _check_out_phase_space_dimension(proc, model, out_phase_space) - - return _differential_cross_section( - proc, model, phase_space_def, in_phase_space, out_phase_space - ) -end - -""" - differential_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, - out_phase_space::AbstractVecOrMat{T}, -) where {T<:Real} - -If the given phase spaces are physical, return differential cross section evaluated at the coordinates. Zero otherwise. -""" -function differential_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, - out_phase_space::AbstractVecOrMat{T}, -) where {T<:Real} - _check_in_phase_space_dimension(proc, model, in_phase_space) - _check_out_phase_space_dimension(proc, model, out_phase_space) - - return _differential_cross_section( - proc, model, phase_space_def, in_phase_space, out_phase_space - ) -end - -""" - differential_cross_section(phase_space_point::PhaseSpacePoint) - -If the given phase spaces are physical, return differential cross section evaluated on a phase space point. Zero otherwise. -""" -function differential_cross_section(phase_space_point::PhaseSpacePoint) - return _differential_cross_section( - phase_space_point.proc, - phase_space_point.model, - phase_space_point.ps_def, - momentum.(phase_space_point.in_particles), - momentum.(phase_space_point.out_particles), - ) -end - -############ -# Total cross sections -############ - -# total cross section on single phase space point -# based on four-momenta -function _total_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVector{T}, -) where {T<:QEDbase.AbstractFourMomentum} - I = 1 / (4 * _incident_flux(proc, model, in_phase_space)) - - return I * _total_probability(proc, model, phase_space_def, in_phase_space) -end - -# total cross section on single phase space point -# based on coordinates -function _total_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVector{T}, -) where {T<:Real} - in_momenta = _generate_incoming_momenta(proc, model, phase_space_def, in_phase_space) - return _total_cross_section(proc, model, phase_space_def, in_momenta) -end - -# total cross section on several phase space points -function _total_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - in_phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractMatrix{T}, -) where {T<:AbstractPhasespaceElement} - res = Vector{eltype(T)}(undef, size(in_phase_space, 2)) - for i in 1:size(in_phase_space, 2) - res[i] = _total_cross_section( - proc, model, in_phase_space_def, view(in_phase_space, :, i) - ) - end - return res -end - -""" - total_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - in_phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, - ) where {T<:QEDbase.AbstractFourMomentum} - -Return the total cross section for a given combination of scattering process and compute model, evaluated at the particle momenta. -""" -function total_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - in_phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, -) where {T<:QEDbase.AbstractFourMomentum} - _check_in_phase_space_dimension(proc, model, in_phase_space) - - return _total_cross_section(proc, model, in_phase_space_def, in_phase_space) -end - -""" - total_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - in_phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, - ) where {T<:Real} - -Return the total cross section for a given combination of scattering process and compute model, evaluated at the coordinates. -""" -function total_cross_section( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - in_phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, -) where {T<:Real} - _check_in_phase_space_dimension(proc, model, in_phase_space) - - return _total_cross_section(proc, model, in_phase_space_def, in_phase_space) -end diff --git a/src/interfaces/process_interface.jl b/src/interfaces/process_interface.jl index 70b5a9a..c2ba09b 100644 --- a/src/interfaces/process_interface.jl +++ b/src/interfaces/process_interface.jl @@ -22,7 +22,11 @@ interface functions need to be implemented for every combination of `CustomProce `CustomModel<:AbstractModelDefinition`, and `CustomPhasespaceDefinition<:AbstractPhasespaceDefinition`. ```Julia - _incident_flux(psp::PhaseSpacePoint{CustomProcess,CustomModel}) + _incident_flux( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + in_phase_space::AbstractVector{T} + ) where {T<:QEDbase.AbstractFourMomentum} _matrix_element(psp::PhaseSpacePoint{CustomProcess,CustomModel}) @@ -67,10 +71,11 @@ This function needs to be given to implement the scattering process interface. function outgoing_particles end """ - _incident_flux(PhaseSpacePoint{PROC,MODEL}) where { - PROC <: AbstractProcessDefinition, - MODEL <: AbstractModelDefinition - } + _incident_flux( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + in_phase_space::AbstractVector{T} + ) where {T<:QEDbase.AbstractFourMomentum} Interface function, which returns the incident flux of the given scattering process for a given incoming phase-space. @@ -87,15 +92,6 @@ Interface function, which returns a tuple of scattering matrix elements for each """ function _matrix_element end -function _matrix_element_square( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - in_phase_space::AbstractVector{T}, - out_phase_space::AbstractVector{T}, -) where {T<:QEDbase.AbstractFourMomentum} - mat_el = _matrix_element(proc, model, in_phase_space, out_phase_space) - return abs2.(mat_el) -end """ _averaging_norm( diff --git a/src/probabilities.jl b/src/probabilities.jl deleted file mode 100644 index cee0181..0000000 --- a/src/probabilities.jl +++ /dev/null @@ -1,290 +0,0 @@ -############ -# scattering probabilities -# -# This file contains implementations of the scattering probability based on the -# process interface with and without input validation and/or phase space -# constraint. -############ - -# differential probability without energy momentum conservation check -# single in phase space points/ single out phase space point -# based on four-momenta -function _unsafe_differential_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVector{T}, - out_phase_space::AbstractVector{T}, -) where {T<:QEDbase.AbstractFourMomentum} - matrix_elements_sq = _matrix_element_square( - proc, model, in_phase_space, out_phase_space - ) - - normalization = _averaging_norm(proc) - - ps_fac = _phase_space_factor( - proc, model, phase_space_def, in_phase_space, out_phase_space - ) - - return normalization * sum(matrix_elements_sq) * ps_fac -end - -# differential probability without energy momentum conservation check -# based on coordinates -function _unsafe_differential_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVector{T}, - out_phase_space::AbstractVector{T}, -) where {T<:Real} - in_momenta, out_momenta = _generate_momenta( - proc, model, phase_space_def, in_phase_space, out_phase_space - ) - return _unsafe_differential_probability( - proc, model, phase_space_def, in_momenta, out_momenta - ) -end - -""" - unsafe_differential_probability(phase_space_point::PhaseSpacePoint) - -Return differential probability evaluated on a phase space point without checking if the given phase space(s) are physical. -""" -function unsafe_differential_probability(phase_space_point::PhaseSpacePoint) - return _unsafe_differential_probability( - phase_space_point.proc, - phase_space_point.model, - phase_space_point.ps_def, - momentum.(phase_space_point.in_particles), - momentum.(phase_space_point.out_particles), - ) -end - -""" - unsafe_differential_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, - out_phase_space::AbstractVecOrMat{T}, - ) where {T<:QEDbase.AbstractFourMomentum} - -Return differential probability evaluated at the four-momenta without checking if the given phase space(s) are physical. -""" -function unsafe_differential_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, - out_phase_space::AbstractVecOrMat{T}, -) where {T<:QEDbase.AbstractFourMomentum} - _check_in_phase_space_dimension(proc, model, in_phase_space) - _check_out_phase_space_dimension(proc, model, out_phase_space) - - return _unsafe_differential_probability( - proc, model, phase_space_def, in_phase_space, out_phase_space - ) -end - -""" - unsafe_differential_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, - out_phase_space::AbstractVecOrMat{T}, - ) where {T<:Real} - -Return differential probability evaluated at the coordinates without checking if the given phase space(s) are physical. -""" -function unsafe_differential_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, - out_phase_space::AbstractVecOrMat{T}, -) where {T<:Real} - _check_in_phase_space_dimension(proc, model, in_phase_space) - _check_out_phase_space_dimension(proc, model, out_phase_space) - - return _unsafe_differential_probability( - proc, model, phase_space_def, in_phase_space, out_phase_space - ) -end - -# differential probability with energy momentum conservation check -# one in phase space point/ one out phase space point -# based on four-momenta -function _differential_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVector{T}, - out_phase_space::AbstractVector{T}, -) where {T<:QEDbase.AbstractFourMomentum} - if !_is_in_phasespace(proc, model, phase_space_def, in_phase_space, out_phase_space) - return zero(eltype(T)) - end - - return _unsafe_differential_probability( - proc, model, phase_space_def, in_phase_space, out_phase_space - ) -end - -# differential probability with energy momentum conservation check -# one in phase space point/ one out phase space point -# based on coordinates -function _differential_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVector{T}, - out_phase_space::AbstractVector{T}, -) where {T<:Real} - in_momenta, out_momenta = _generate_momenta( - proc, model, phase_space_def, in_phase_space, out_phase_space - ) - return _differential_probability(proc, model, phase_space_def, in_momenta, out_momenta) -end - -""" - differential_probability(phase_space_point::PhaseSpacePoint) - -If the given phase spaces are physical, return differential probability evaluated on a phase space point. Zero otherwise. -""" -function differential_probability(phase_space_point::PhaseSpacePoint) - return differential_probability( - phase_space_point.proc, - phase_space_point.model, - phase_space_point.ps_def, - momentum.(phase_space_point.in_particles), - momentum.(phase_space_point.out_particles), - ) -end - -""" - differential_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, - out_phase_space::AbstractVecOrMat{T}, - ) where {T<:QEDbase.AbstractFourMomentum} - -If the given phase spaces are physical, return differential probability evaluated at the four-momenta. Zero otherwise. -""" -function differential_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, - out_phase_space::AbstractVecOrMat{T}, -) where {T<:QEDbase.AbstractFourMomentum} - _check_in_phase_space_dimension(proc, model, in_phase_space) - _check_out_phase_space_dimension(proc, model, out_phase_space) - - return _differential_probability( - proc, model, phase_space_def, in_phase_space, out_phase_space - ) -end - -""" - differential_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, - out_phase_space::AbstractVecOrMat{T}, -) where {T<:Real} - -If the given phase spaces are physical, return differential probability evaluated at the coordinates. Zero otherwise. -""" -function differential_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, - out_phase_space::AbstractVecOrMat{T}, -) where {T<:Real} - _check_in_phase_space_dimension(proc, model, in_phase_space) - _check_out_phase_space_dimension(proc, model, out_phase_space) - - return _differential_probability( - proc, model, phase_space_def, in_phase_space, out_phase_space - ) -end - -########### -# Total probability -########### - -# total probability on a phase space point -# based on coordinates -function _total_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVector{T}, -) where {T<:Real} - in_momenta = _generate_incoming_momenta(proc, model, phase_space_def, in_phase_space) - return _total_probability(proc, model, phase_space_def, in_momenta) -end - -# total probability on several phase space points -function _total_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractMatrix{T}, -) where {T<:AbstractPhasespaceElement} - res = Vector{eltype(T)}(undef, size(in_phase_space, 2)) - for i in 1:size(in_phase_space, 2) - res[i] = _total_probability( - proc, model, phase_space_def, view(in_phase_space, :, i) - ) - end - return res -end - -""" - total_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractMatrix{T}, - ) where {T<:QEDbase.AbstractFourMomentum} - -Return the total probability of a given model and process combination, evaluated at the particle momenta. -""" -function total_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, -) where {T<:QEDbase.AbstractFourMomentum} - _check_in_phase_space_dimension(proc, model, in_phase_space) - - return _total_probability(proc, model, phase_space_def, in_phase_space) -end - -""" - total_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractMatrix{T}, - ) where {T<:Real} - -Return the total probability of a given model and process combination, evaluated at the coordinates. -""" -function total_probability( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - phase_space_def::AbstractPhasespaceDefinition, - in_phase_space::AbstractVecOrMat{T}, -) where {T<:Real} - _check_in_phase_space_dimension(proc, model, in_phase_space) - - return _total_probability(proc, model, phase_space_def, in_phase_space) -end diff --git a/test/cross_sections.jl b/test/cross_sections.jl index 79d832d..40f9921 100644 --- a/test/cross_sections.jl +++ b/test/cross_sections.jl @@ -26,13 +26,11 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() p_in_unphys_invalid = TestImplementation._rand_in_momenta_failing(RNG, N_INCOMING + 1) p_out_phys = TestImplementation._rand_momenta(RNG, N_OUTGOING) - p_out_phys_invalid = TestImplementation._rand_momenta(RNG, N_OUTGOING + 1) p_out_unphys = TestImplementation._rand_out_momenta_failing(RNG, N_OUTGOING) - p_out_unphys_invalid = TestImplementation._rand_out_momenta_failing(RNG, N_OUTGOING + 1) - p_in_all = (p_in_phys, p_in_unphys, p_in_phys_invalid, p_in_unphys_invalid) + p_in_all = (p_in_phys, p_in_unphys) - p_out_all = (p_out_phys, p_out_phys_invalid, p_out_unphys, p_out_unphys_invalid) + p_out_all = (p_out_phys, p_out_unphys) # all combinations p_combs = Iterators.product(p_in_all, p_out_all) @@ -41,6 +39,8 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() p_out_all_valid = (p_out_phys, p_out_unphys) + p_in_all_valid = (p_in_phys, p_in_unphys ) + # all valid combinations p_combs_valid = Iterators.product(p_in_all_valid, p_out_all_valid) @@ -50,229 +50,90 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() p_combs_phys = Iterators.product(p_in_all_phys, p_out_all_phys) @testset "cross section" begin - @testset "unsafe" begin - @testset "compute" begin - for (P_IN, P_OUT) in p_combs_phys - diffCS = unsafe_differential_cross_section( - TESTPROC, TESTMODEL, TESTPSDEF, P_IN, P_OUT - ) - groundtruth = TestImplementation._groundtruth_unsafe_diffCS( - TESTPROC, P_IN, P_OUT - ) - @test isapprox(diffCS, groundtruth, atol=ATOL, rtol=RTOL) - end - end - - @testset "compute on phase space points" begin + @testset "unsafe compute" begin + for (P_IN, P_OUT) in p_combs_phys PS_POINT = generate_phase_space( - TESTPROC, TESTMODEL, TESTPSDEF, p_in_phys, p_out_phys + TESTPROC, TESTMODEL, TESTPSDEF,P_IN, P_OUT ) diffCS_on_psp = unsafe_differential_cross_section(PS_POINT) groundtruth = TestImplementation._groundtruth_unsafe_diffCS( - TESTPROC, p_in_phys, p_out_phys + TESTPROC, P_IN, P_OUT ) @test isapprox(diffCS_on_psp, groundtruth, atol=ATOL, rtol=RTOL) end + end - @testset "invalid input" begin - for (P_IN, P_OUT) in p_combs - - # filter out all valid combinations - if !((P_IN, P_OUT) in p_combs_valid) - @test_throws DimensionMismatch unsafe_differential_cross_section( - TESTPROC, TESTMODEL, TESTPSDEF, P_IN, P_OUT - ) - - COORDS_IN = TestImplementation.flat_components(P_IN) - COORDS_OUT = TestImplementation.flat_components(P_OUT) - @test_throws DimensionMismatch unsafe_differential_cross_section( - TESTPROC, TESTMODEL, TESTPSDEF, COORDS_IN, COORDS_OUT - ) - end - end + @testset "safe compute" begin + for (P_IN, P_OUT) in p_combs_valid + PS_POINT = generate_phase_space( + TESTPROC, TESTMODEL, TESTPSDEF,P_IN, P_OUT + ) + diffCS_on_psp = differential_cross_section(PS_POINT) + groundtruth = TestImplementation._groundtruth_safe_diffCS( + TESTPROC, P_IN, P_OUT + ) + @test isapprox(diffCS_on_psp, groundtruth, atol=ATOL, rtol=RTOL) end + end - @testset "safe" begin - @testset "compute" begin - for (P_IN, P_OUT) in p_combs_valid - diffCS_on_moms = differential_cross_section( - TESTPROC, TESTMODEL, TESTPSDEF, P_IN, P_OUT - ) - - COORDS_IN = TestImplementation.flat_components(P_IN) - COORDS_OUT = TestImplementation.flat_components(P_OUT) - diffCS_on_coords = differential_cross_section( - TESTPROC, TESTMODEL, TESTPSDEF, COORDS_IN, COORDS_OUT - ) - groundtruth = TestImplementation._groundtruth_safe_diffCS( - TESTPROC, P_IN, P_OUT - ) - @test isapprox(diffCS_on_moms, groundtruth, atol=ATOL, rtol=RTOL) - @test isapprox(diffCS_on_coords, groundtruth, atol=ATOL, rtol=RTOL) - end - end - - @testset "compute on phase space points" begin - PS_POINT = generate_phase_space( - TESTPROC, TESTMODEL, TESTPSDEF, p_in_phys, p_out_phys + @testset "total cross section" begin + @testset "compute" begin + for P_IN in (p_in_phys,) + groundtruth = TestImplementation._groundtruth_total_cross_section( + P_IN ) - diffCS_on_psp = differential_cross_section(PS_POINT) - groundtruth = TestImplementation._groundtruth_safe_diffCS( - TESTPROC, p_in_phys, p_out_phys + totCS_on_moms = total_cross_section( + TESTPROC, TESTMODEL, TESTPSDEF, P_IN ) - @test isapprox(diffCS_on_psp, groundtruth, atol=ATOL, rtol=RTOL) - end - - @testset "invalid input" begin - for (P_IN, P_OUT) in p_combs - # filter out all valid combinations - if !((P_IN, P_OUT) in p_combs_valid) - @test_throws DimensionMismatch differential_cross_section( - TESTPROC, TESTMODEL, TESTPSDEF, P_IN, P_OUT - ) + COORDS_IN = TestImplementation.flat_components(P_IN) + totCS_on_coords = total_cross_section( + TESTPROC, TESTMODEL, TESTPSDEF, COORDS_IN + ) - COORDS_IN = TestImplementation.flat_components(P_IN) - COORDS_OUT = TestImplementation.flat_components(P_OUT) - @test_throws DimensionMismatch differential_cross_section( - TESTPROC, TESTMODEL, TESTPSDEF, COORDS_IN, COORDS_OUT - ) - end - end + @test isapprox(totCS_on_moms, groundtruth, atol=ATOL, rtol=RTOL) + @test isapprox(totCS_on_coords, groundtruth, atol=ATOL, rtol=RTOL) end end + @testset "invalid input" begin + for P_IN in (p_in_phys_invalid,) + @test_throws DimensionMismatch total_cross_section( + TESTPROC, TESTMODEL, TESTPSDEF, P_IN + ) - @testset "total cross section" begin - @testset "compute" begin - for P_IN in (p_in_phys,) - groundtruth = TestImplementation._groundtruth_total_cross_section( - P_IN - ) - totCS_on_moms = total_cross_section( - TESTPROC, TESTMODEL, TESTPSDEF, P_IN - ) - - COORDS_IN = TestImplementation.flat_components(P_IN) - totCS_on_coords = total_cross_section( - TESTPROC, TESTMODEL, TESTPSDEF, COORDS_IN - ) - - @test isapprox(totCS_on_moms, groundtruth, atol=ATOL, rtol=RTOL) - @test isapprox(totCS_on_coords, groundtruth, atol=ATOL, rtol=RTOL) - end - end - @testset "invalid input" begin - for P_IN in (p_in_phys_invalid,) - @test_throws DimensionMismatch total_cross_section( - TESTPROC, TESTMODEL, TESTPSDEF, P_IN - ) - - COORDS_IN = TestImplementation.flat_components(P_IN) - @test_throws DimensionMismatch total_cross_section( - TESTPROC, TESTMODEL, TESTPSDEF, COORDS_IN - ) - end + COORDS_IN = TestImplementation.flat_components(P_IN) + @test_throws DimensionMismatch total_cross_section( + TESTPROC, TESTMODEL, TESTPSDEF, COORDS_IN + ) end end end @testset "differential probability" begin - @testset "unsafe" begin - @testset "compute" begin - for (P_IN, P_OUT) in p_combs_phys - prob_on_moms = unsafe_differential_probability( - TESTPROC, TESTMODEL, TESTPSDEF, P_IN, P_OUT - ) - COORDS_IN = TestImplementation.flat_components(P_IN) - COORDS_OUT = TestImplementation.flat_components(P_OUT) - prob_on_coords = unsafe_differential_probability( - TESTPROC, TESTMODEL, TESTPSDEF, COORDS_IN, COORDS_OUT - ) - groundtruth = TestImplementation._groundtruth_unsafe_probability( - TESTPROC, P_IN, P_OUT - ) - @test isapprox(prob_on_moms, groundtruth, atol=ATOL, rtol=RTOL) - @test isapprox(prob_on_coords, groundtruth, atol=ATOL, rtol=RTOL) - end - end - - @testset "compute on phase space points" begin + @testset "unsafe compute" begin + for (P_IN, P_OUT) in p_combs_phys PS_POINT = generate_phase_space( - TESTPROC, TESTMODEL, TESTPSDEF, p_in_phys, p_out_phys + TESTPROC, TESTMODEL, TESTPSDEF, P_IN, P_OUT ) prop_on_psp = unsafe_differential_probability(PS_POINT) groundtruth = TestImplementation._groundtruth_unsafe_probability( - TESTPROC, p_in_phys, p_out_phys + TESTPROC, P_IN, P_OUT ) @test isapprox(prop_on_psp, groundtruth, atol=ATOL, rtol=RTOL) end - - @testset "invalid input" begin - for (P_IN, P_OUT) in p_combs - - # filter out all valid combinations - if !((P_IN, P_OUT) in p_combs_valid) - @test_throws DimensionMismatch unsafe_differential_probability( - TESTPROC, TESTMODEL, TESTPSDEF, P_IN, P_OUT - ) - - COORDS_IN = TestImplementation.flat_components(P_IN) - COORDS_OUT = TestImplementation.flat_components(P_OUT) - @test_throws DimensionMismatch unsafe_differential_probability( - TESTPROC, TESTMODEL, TESTPSDEF, COORDS_IN, COORDS_OUT - ) - end - end - end end - @testset "safe" begin - @testset "compute" begin - for (P_IN, P_OUT) in p_combs_valid - prob_on_moms = differential_probability( - TESTPROC, TESTMODEL, TESTPSDEF, P_IN, P_OUT - ) - - COORDS_IN = TestImplementation.flat_components(P_IN) - COORDS_OUT = TestImplementation.flat_components(P_OUT) - prob_on_coords = differential_probability( - TESTPROC, TESTMODEL, TESTPSDEF, COORDS_IN, COORDS_OUT - ) - groundtruth = TestImplementation._groundtruth_safe_probability( - TESTPROC, P_IN, P_OUT - ) - @test isapprox(prob_on_moms, groundtruth, atol=ATOL, rtol=RTOL) - @test isapprox(prob_on_coords, groundtruth, atol=ATOL, rtol=RTOL) - end - end - @testset "compute on phase space points" begin + @testset "safe compute" begin + for (P_IN, P_OUT) in p_combs_valid PS_POINT = generate_phase_space( - TESTPROC, TESTMODEL, TESTPSDEF, p_in_phys, p_out_phys + TESTPROC, TESTMODEL, TESTPSDEF, P_IN, P_OUT ) prop_on_psp = differential_probability(PS_POINT) groundtruth = TestImplementation._groundtruth_safe_probability( - TESTPROC, p_in_phys, p_out_phys + TESTPROC, P_IN, P_OUT ) @test isapprox(prop_on_psp, groundtruth, atol=ATOL, rtol=RTOL) end - - @testset "invalid input" begin - for (P_IN, P_OUT) in p_combs - - # filter out all valid combinations - if !((P_IN, P_OUT) in p_combs_valid) - @test_throws DimensionMismatch differential_probability( - TESTPROC, TESTMODEL, TESTPSDEF, P_IN, P_OUT - ) - - COORDS_IN = TestImplementation.flat_components(P_IN) - COORDS_OUT = TestImplementation.flat_components(P_OUT) - @test_throws DimensionMismatch differential_probability( - TESTPROC, TESTMODEL, TESTPSDEF, COORDS_IN, COORDS_OUT - ) - end - end - end end @testset "total probability" begin @@ -295,7 +156,7 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() end end @testset "invalid input" begin - for P_IN in (p_in_phys_invalid) + for P_IN in (p_in_phys_invalid,) @test_throws DimensionMismatch total_probability( TESTPROC, TESTMODEL, TESTPSDEF, P_IN ) diff --git a/test/interfaces/process_interface.jl b/test/interfaces/process_interface.jl index 445ceb3..960dc6d 100644 --- a/test/interfaces/process_interface.jl +++ b/test/interfaces/process_interface.jl @@ -72,7 +72,7 @@ include("../test_implementation/TestImplementation.jl") end @testset "incident flux" begin - test_incident_flux = QEDprocesses._incident_flux(PSP) + test_incident_flux = QEDprocesses._incident_flux(TESTPROC,TESTMODEL,IN_PS) groundtruth = TestImplementation._groundtruth_incident_flux(IN_PS) @test isapprox(test_incident_flux, groundtruth, atol=ATOL, rtol=RTOL) end diff --git a/test/test_implementation/test_process.jl b/test/test_implementation/test_process.jl index 30f9781..8dfc209 100644 --- a/test/test_implementation/test_process.jl +++ b/test/test_implementation/test_process.jl @@ -71,8 +71,9 @@ struct TestPhasespaceDef_FAIL <: AbstractPhasespaceDefinition end # dummy implementation of the process interface -function QEDprocesses._incident_flux(psp::PhaseSpacePoint{<:TestProcess,TestModel}) - in_ps = momentum.(psp.in_particles) +function QEDprocesses._incident_flux( + ::TestProcess, ::TestModel, in_ps::AbstractVector{T} +) where {T<:QEDbase.AbstractFourMomentum} return _groundtruth_incident_flux(in_ps) end @@ -118,7 +119,8 @@ function QEDprocesses._generate_outgoing_momenta( return _groundtruth_generate_momenta(out_phase_space) end -function QEDprocesses._total_probability(psp::PhaseSpacePoint{<:TestProcess,TestModel}) - in_ps = momentum.(psp.in_particles) +function QEDprocesses._total_probability( + proc::TestProcess, model::TestModel, ps_def::TestPhasespaceDef, in_ps::AbstractVector{T} +) where {T<:QEDbase.AbstractFourMomentum} return _groundtruth_total_probability(in_ps) end From dff0d2007c5c7e4b0a6f6046537b0f39eb3e5371 Mon Sep 17 00:00:00 2001 From: Uwe Hernandez Acosta Date: Thu, 16 May 2024 23:52:30 +0200 Subject: [PATCH 07/22] updated perturbative Compton to PSPs --- src/cross_section/internal.jl | 2 +- .../perturbative/cross_section.jl | 35 +++++++------------ .../one_photon_compton/perturbative.jl | 35 +++++++++++++------ test/runtests.jl | 2 +- 4 files changed, 39 insertions(+), 35 deletions(-) diff --git a/src/cross_section/internal.jl b/src/cross_section/internal.jl index d0e8f96..a8bf22d 100644 --- a/src/cross_section/internal.jl +++ b/src/cross_section/internal.jl @@ -76,7 +76,7 @@ function _unsafe_differential_cross_section( out_phase_space::AbstractVector{T}, ) where {T<:QEDbase.AbstractFourMomentum} psp = generate_phase_space(proc,model,phase_space_def,in_phase_space,out_phase_space) - return _unsafe_differential_cross_section(phase_space_point) + return unsafe_differential_cross_section(psp) end # convenience function for internal use only diff --git a/src/processes/one_photon_compton/perturbative/cross_section.jl b/src/processes/one_photon_compton/perturbative/cross_section.jl index 4330433..2b4bd0f 100644 --- a/src/processes/one_photon_compton/perturbative/cross_section.jl +++ b/src/processes/one_photon_compton/perturbative/cross_section.jl @@ -9,13 +9,10 @@ function _incident_flux( return @inbounds in_ps[1] * in_ps[2] end -function _matrix_element( - proc::Compton, - model::PerturbativeQED, - in_ps::AbstractVector{T}, - out_ps::AbstractVector{T}, -) where {T<:QEDbase.AbstractFourMomentum} - return _pert_compton_matrix_element(proc, in_ps, out_ps) +function _matrix_element(psp::PhaseSpacePoint{<:Compton,PerturbativeQED}) + in_ps = momentum.(psp.in_particles) + out_ps = momentum.(psp.out_particles) + return _pert_compton_matrix_element(psp.proc, in_ps, out_ps) end """ @@ -40,24 +37,16 @@ function _all_onshell( isapprox(sq_out_moms, SVector(sq_out_masses)) end -function _is_in_phasespace( - proc::Compton, - model::PerturbativeQED, - in_ps_def::AbstractPhasespaceDefinition, - in_ps::AbstractVector{T}, - out_ps::AbstractVector{T}, -) where {T<:QEDbase.AbstractFourMomentum} - return (isapprox(sum(in_ps), sum(out_ps))) ? _all_onshell(proc, in_ps, out_ps) : false +function _is_in_phasespace(psp::PhaseSpacePoint{<:Compton,PerturbativeQED}) + in_ps = momentum.(psp.in_particles) + out_ps = momentum.(psp.out_particles) + return (isapprox(sum(in_ps), sum(out_ps))) ? _all_onshell(psp.proc, in_ps, out_ps) : false end -@inline function _phase_space_factor( - proc::Compton, - model::PerturbativeQED, - in_ps_def::AbstractPhasespaceDefinition, - in_ps::AbstractVector{T}, - out_ps::AbstractVector{T}, -) where {T<:QEDbase.AbstractFourMomentum} - return _pert_compton_ps_fac(in_ps_def, in_ps[2], out_ps[2]) +@inline function _phase_space_factor(psp::PhaseSpacePoint{<:Compton,PerturbativeQED}) + in_ps = momentum.(psp.in_particles) + out_ps = momentum.(psp.out_particles) + return _pert_compton_ps_fac(psp.ps_def, in_ps[2], out_ps[2]) end ####### diff --git a/test/processes/one_photon_compton/perturbative.jl b/test/processes/one_photon_compton/perturbative.jl index 9d28869..9258269 100644 --- a/test/processes/one_photon_compton/perturbative.jl +++ b/test/processes/one_photon_compton/perturbative.jl @@ -48,12 +48,17 @@ end Iterators.product(COS_THETAS, PHIS) IN_COORDS = [omega] OUT_COORDS = [cos_theta, phi] + IN_PS, OUT_PS = QEDprocesses._generate_momenta( + PROC, MODEL, PS_DEF, IN_COORDS, OUT_COORDS + ) + groundtruth = _groundtruth_pert_compton_diffCS_spinsum_polsum( omega, cos_theta, 1.0 ) - test_val = unsafe_differential_cross_section( - PROC, MODEL, PS_DEF, IN_COORDS, OUT_COORDS - ) + + PSP = generate_phase_space(PROC,MODEL,PS_DEF,IN_PS,OUT_PS) + test_val = unsafe_differential_cross_section(PSP) + @test isapprox(test_val, groundtruth, atol=ATOL, rtol=RTOL) end end @@ -65,12 +70,17 @@ end Iterators.product(COS_THETAS, PHIS) IN_COORDS = [omega] OUT_COORDS = [cos_theta, phi] + IN_PS, OUT_PS = QEDprocesses._generate_momenta( + PROC, MODEL, PS_DEF, IN_COORDS, OUT_COORDS + ) + groundtruth = _groundtruth_pert_compton_diffCS_spinsum_xpol( omega, cos_theta, phi, 1.0 ) - test_val = unsafe_differential_cross_section( - PROC, MODEL, PS_DEF, IN_COORDS, OUT_COORDS - ) + + PSP = generate_phase_space(PROC,MODEL,PS_DEF,IN_PS,OUT_PS) + test_val = unsafe_differential_cross_section(PSP) + @test isapprox(test_val, groundtruth, atol=ATOL, rtol=RTOL) end end @@ -82,12 +92,17 @@ end Iterators.product(COS_THETAS, PHIS) IN_COORDS = [omega] OUT_COORDS = [cos_theta, phi] + IN_PS, OUT_PS = QEDprocesses._generate_momenta( + PROC, MODEL, PS_DEF, IN_COORDS, OUT_COORDS + ) + groundtruth = _groundtruth_pert_compton_diffCS_spinsum_ypol( omega, cos_theta, phi, 1.0 ) - test_val = unsafe_differential_cross_section( - PROC, MODEL, PS_DEF, IN_COORDS, OUT_COORDS - ) + + PSP = generate_phase_space(PROC,MODEL,PS_DEF,IN_PS,OUT_PS) + test_val = unsafe_differential_cross_section(PSP) + @test isapprox(test_val, groundtruth, atol=ATOL, rtol=RTOL) end end @@ -98,7 +113,7 @@ end # Klein-Nishina: total cross section function klein_nishina_total_cross_section(in_ps) function func(x) - return unsafe_differential_cross_section( + return QEDprocesses._unsafe_differential_cross_section( Compton(), PerturbativeQED(), PS_DEF, in_ps, [x, 0.0] ) end diff --git a/test/runtests.jl b/test/runtests.jl index 76c49ae..9ecc6ee 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -19,7 +19,7 @@ begin include("propagators.jl") end @time @safetestset "cross section & probability" begin - #include("cross_sections.jl") + include("cross_sections.jl") end @time @safetestset "phase spaces" begin From 62c430d5349692bc719a57b7d90b0232947d0a83 Mon Sep 17 00:00:00 2001 From: Uwe Hernandez Acosta Date: Thu, 16 May 2024 23:52:54 +0200 Subject: [PATCH 08/22] formatting --- src/cross_section/diff_cross_section.jl | 13 ++------- src/cross_section/diff_probability.jl | 4 +-- src/cross_section/internal.jl | 16 ++++++++--- src/interfaces/process_interface.jl | 1 - .../perturbative/cross_section.jl | 6 +++- test/cross_sections.jl | 28 ++++++++----------- test/interfaces/process_interface.jl | 2 +- .../one_photon_compton/perturbative.jl | 6 ++-- 8 files changed, 35 insertions(+), 41 deletions(-) diff --git a/src/cross_section/diff_cross_section.jl b/src/cross_section/diff_cross_section.jl index 74bcb8d..24f09bd 100644 --- a/src/cross_section/diff_cross_section.jl +++ b/src/cross_section/diff_cross_section.jl @@ -5,13 +5,8 @@ # cross sections based on the scattering process interface ######################## - function _incident_flux(psp::PhaseSpacePoint) - return _incident_flux( - psp.proc, - psp.model, - momentum.(psp.in_particles) - ) + return _incident_flux(psp.proc, psp.model, momentum.(psp.in_particles)) end """ @@ -25,7 +20,6 @@ function unsafe_differential_cross_section(phase_space_point::PhaseSpacePoint) return I * unsafe_differential_probability(phase_space_point) end - """ differential_cross_section(phase_space_point::PhaseSpacePoint) @@ -33,11 +27,8 @@ If the given phase spaces are physical, return differential cross section evalua """ function differential_cross_section(phase_space_point::PhaseSpacePoint) if !_is_in_phasespace(phase_space_point) - return zero(eltype(momentum(phase_space_point,Incoming(),1))) + return zero(eltype(momentum(phase_space_point, Incoming(), 1))) end return unsafe_differential_cross_section(phase_space_point) end - - - diff --git a/src/cross_section/diff_probability.jl b/src/cross_section/diff_probability.jl index 423f865..9d74dc5 100644 --- a/src/cross_section/diff_probability.jl +++ b/src/cross_section/diff_probability.jl @@ -35,10 +35,8 @@ If the given phase spaces are physical, return differential probability evaluate """ function differential_probability(phase_space_point::PhaseSpacePoint) if !_is_in_phasespace(phase_space_point) - return zero(eltype(momentum(phase_space_point,Incoming(),1))) + return zero(eltype(momentum(phase_space_point, Incoming(), 1))) end return unsafe_differential_probability(phase_space_point) end - - diff --git a/src/cross_section/internal.jl b/src/cross_section/internal.jl index a8bf22d..2dc2a1b 100644 --- a/src/cross_section/internal.jl +++ b/src/cross_section/internal.jl @@ -10,7 +10,9 @@ function _unsafe_differential_probability( in_phase_space::AbstractVector{T}, out_phase_space::AbstractVector{T}, ) where {T<:QEDbase.AbstractFourMomentum} - psp = generate_phase_space(proc,model,phase_space_def,in_phase_space,out_phase_space) + psp = generate_phase_space( + proc, model, phase_space_def, in_phase_space, out_phase_space + ) return unsafe_differential_probability(psp) end @@ -43,7 +45,9 @@ function _differential_probability( in_phase_space::AbstractVector{T}, out_phase_space::AbstractVector{T}, ) where {T<:QEDbase.AbstractFourMomentum} - psp = generate_phase_space(proc,model,phase_space_def,in_phase_space,out_phase_space) + psp = generate_phase_space( + proc, model, phase_space_def, in_phase_space, out_phase_space + ) return differential_probability(psp) end @@ -75,7 +79,9 @@ function _unsafe_differential_cross_section( in_phase_space::AbstractVector{T}, out_phase_space::AbstractVector{T}, ) where {T<:QEDbase.AbstractFourMomentum} - psp = generate_phase_space(proc,model,phase_space_def,in_phase_space,out_phase_space) + psp = generate_phase_space( + proc, model, phase_space_def, in_phase_space, out_phase_space + ) return unsafe_differential_cross_section(psp) end @@ -108,7 +114,9 @@ function _differential_cross_section( in_phase_space::AbstractVector{T}, out_phase_space::AbstractVector{T}, ) where {T<:QEDbase.AbstractFourMomentum} - psp = generate_phase_space(proc,model,phase_space_def,in_phase_space,out_phase_space) + psp = generate_phase_space( + proc, model, phase_space_def, in_phase_space, out_phase_space + ) return differential_cross_section(psp) end diff --git a/src/interfaces/process_interface.jl b/src/interfaces/process_interface.jl index c2ba09b..3c6c0e4 100644 --- a/src/interfaces/process_interface.jl +++ b/src/interfaces/process_interface.jl @@ -92,7 +92,6 @@ Interface function, which returns a tuple of scattering matrix elements for each """ function _matrix_element end - """ _averaging_norm( proc::AbstractProcessDefinition diff --git a/src/processes/one_photon_compton/perturbative/cross_section.jl b/src/processes/one_photon_compton/perturbative/cross_section.jl index 2b4bd0f..c4a5891 100644 --- a/src/processes/one_photon_compton/perturbative/cross_section.jl +++ b/src/processes/one_photon_compton/perturbative/cross_section.jl @@ -40,7 +40,11 @@ end function _is_in_phasespace(psp::PhaseSpacePoint{<:Compton,PerturbativeQED}) in_ps = momentum.(psp.in_particles) out_ps = momentum.(psp.out_particles) - return (isapprox(sum(in_ps), sum(out_ps))) ? _all_onshell(psp.proc, in_ps, out_ps) : false + return if (isapprox(sum(in_ps), sum(out_ps))) + _all_onshell(psp.proc, in_ps, out_ps) + else + false + end end @inline function _phase_space_factor(psp::PhaseSpacePoint{<:Compton,PerturbativeQED}) diff --git a/test/cross_sections.jl b/test/cross_sections.jl index 40f9921..c99cefc 100644 --- a/test/cross_sections.jl +++ b/test/cross_sections.jl @@ -28,9 +28,9 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() p_out_phys = TestImplementation._rand_momenta(RNG, N_OUTGOING) p_out_unphys = TestImplementation._rand_out_momenta_failing(RNG, N_OUTGOING) - p_in_all = (p_in_phys, p_in_unphys) + p_in_all = (p_in_phys, p_in_unphys) - p_out_all = (p_out_phys, p_out_unphys) + p_out_all = (p_out_phys, p_out_unphys) # all combinations p_combs = Iterators.product(p_in_all, p_out_all) @@ -39,7 +39,7 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() p_out_all_valid = (p_out_phys, p_out_unphys) - p_in_all_valid = (p_in_phys, p_in_unphys ) + p_in_all_valid = (p_in_phys, p_in_unphys) # all valid combinations p_combs_valid = Iterators.product(p_in_all_valid, p_out_all_valid) @@ -52,12 +52,10 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() @testset "cross section" begin @testset "unsafe compute" begin for (P_IN, P_OUT) in p_combs_phys - PS_POINT = generate_phase_space( - TESTPROC, TESTMODEL, TESTPSDEF,P_IN, P_OUT - ) + PS_POINT = generate_phase_space(TESTPROC, TESTMODEL, TESTPSDEF, P_IN, P_OUT) diffCS_on_psp = unsafe_differential_cross_section(PS_POINT) groundtruth = TestImplementation._groundtruth_unsafe_diffCS( - TESTPROC, P_IN, P_OUT + TESTPROC, P_IN, P_OUT ) @test isapprox(diffCS_on_psp, groundtruth, atol=ATOL, rtol=RTOL) end @@ -65,12 +63,10 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() @testset "safe compute" begin for (P_IN, P_OUT) in p_combs_valid - PS_POINT = generate_phase_space( - TESTPROC, TESTMODEL, TESTPSDEF,P_IN, P_OUT - ) + PS_POINT = generate_phase_space(TESTPROC, TESTMODEL, TESTPSDEF, P_IN, P_OUT) diffCS_on_psp = differential_cross_section(PS_POINT) groundtruth = TestImplementation._groundtruth_safe_diffCS( - TESTPROC, P_IN, P_OUT + TESTPROC, P_IN, P_OUT ) @test isapprox(diffCS_on_psp, groundtruth, atol=ATOL, rtol=RTOL) end @@ -79,9 +75,7 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() @testset "total cross section" begin @testset "compute" begin for P_IN in (p_in_phys,) - groundtruth = TestImplementation._groundtruth_total_cross_section( - P_IN - ) + groundtruth = TestImplementation._groundtruth_total_cross_section(P_IN) totCS_on_moms = total_cross_section( TESTPROC, TESTMODEL, TESTPSDEF, P_IN ) @@ -113,11 +107,11 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() @testset "unsafe compute" begin for (P_IN, P_OUT) in p_combs_phys PS_POINT = generate_phase_space( - TESTPROC, TESTMODEL, TESTPSDEF, P_IN, P_OUT + TESTPROC, TESTMODEL, TESTPSDEF, P_IN, P_OUT ) prop_on_psp = unsafe_differential_probability(PS_POINT) groundtruth = TestImplementation._groundtruth_unsafe_probability( - TESTPROC, P_IN, P_OUT + TESTPROC, P_IN, P_OUT ) @test isapprox(prop_on_psp, groundtruth, atol=ATOL, rtol=RTOL) end @@ -126,7 +120,7 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() @testset "safe compute" begin for (P_IN, P_OUT) in p_combs_valid PS_POINT = generate_phase_space( - TESTPROC, TESTMODEL, TESTPSDEF, P_IN, P_OUT + TESTPROC, TESTMODEL, TESTPSDEF, P_IN, P_OUT ) prop_on_psp = differential_probability(PS_POINT) groundtruth = TestImplementation._groundtruth_safe_probability( diff --git a/test/interfaces/process_interface.jl b/test/interfaces/process_interface.jl index 960dc6d..d081930 100644 --- a/test/interfaces/process_interface.jl +++ b/test/interfaces/process_interface.jl @@ -72,7 +72,7 @@ include("../test_implementation/TestImplementation.jl") end @testset "incident flux" begin - test_incident_flux = QEDprocesses._incident_flux(TESTPROC,TESTMODEL,IN_PS) + test_incident_flux = QEDprocesses._incident_flux(TESTPROC, TESTMODEL, IN_PS) groundtruth = TestImplementation._groundtruth_incident_flux(IN_PS) @test isapprox(test_incident_flux, groundtruth, atol=ATOL, rtol=RTOL) end diff --git a/test/processes/one_photon_compton/perturbative.jl b/test/processes/one_photon_compton/perturbative.jl index 9258269..cf5ae27 100644 --- a/test/processes/one_photon_compton/perturbative.jl +++ b/test/processes/one_photon_compton/perturbative.jl @@ -56,7 +56,7 @@ end omega, cos_theta, 1.0 ) - PSP = generate_phase_space(PROC,MODEL,PS_DEF,IN_PS,OUT_PS) + PSP = generate_phase_space(PROC, MODEL, PS_DEF, IN_PS, OUT_PS) test_val = unsafe_differential_cross_section(PSP) @test isapprox(test_val, groundtruth, atol=ATOL, rtol=RTOL) @@ -78,7 +78,7 @@ end omega, cos_theta, phi, 1.0 ) - PSP = generate_phase_space(PROC,MODEL,PS_DEF,IN_PS,OUT_PS) + PSP = generate_phase_space(PROC, MODEL, PS_DEF, IN_PS, OUT_PS) test_val = unsafe_differential_cross_section(PSP) @test isapprox(test_val, groundtruth, atol=ATOL, rtol=RTOL) @@ -100,7 +100,7 @@ end omega, cos_theta, phi, 1.0 ) - PSP = generate_phase_space(PROC,MODEL,PS_DEF,IN_PS,OUT_PS) + PSP = generate_phase_space(PROC, MODEL, PS_DEF, IN_PS, OUT_PS) test_val = unsafe_differential_cross_section(PSP) @test isapprox(test_val, groundtruth, atol=ATOL, rtol=RTOL) From 63c3ec6ddff1b7f387f4897c2c7a47aacc52e86f Mon Sep 17 00:00:00 2001 From: Uwe Hernandez Acosta Date: Fri, 17 May 2024 00:19:48 +0200 Subject: [PATCH 09/22] cleanup tests for cross sections and probabilities; formatting --- test/cross_sections.jl | 185 +++++++++++++++++------------------------ 1 file changed, 76 insertions(+), 109 deletions(-) diff --git a/test/cross_sections.jl b/test/cross_sections.jl index c99cefc..f7e6215 100644 --- a/test/cross_sections.jl +++ b/test/cross_sections.jl @@ -23,7 +23,6 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() p_in_phys = TestImplementation._rand_momenta(RNG, N_INCOMING) p_in_phys_invalid = TestImplementation._rand_momenta(RNG, N_INCOMING + 1) p_in_unphys = TestImplementation._rand_in_momenta_failing(RNG, N_INCOMING) - p_in_unphys_invalid = TestImplementation._rand_in_momenta_failing(RNG, N_INCOMING + 1) p_out_phys = TestImplementation._rand_momenta(RNG, N_OUTGOING) p_out_unphys = TestImplementation._rand_out_momenta_failing(RNG, N_OUTGOING) @@ -35,133 +34,101 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() # all combinations p_combs = Iterators.product(p_in_all, p_out_all) - p_in_all_valid = (p_in_phys, p_in_unphys) - - p_out_all_valid = (p_out_phys, p_out_unphys) - - p_in_all_valid = (p_in_phys, p_in_unphys) - - # all valid combinations - p_combs_valid = Iterators.product(p_in_all_valid, p_out_all_valid) - - p_in_all_phys = (p_in_phys,) - p_out_all_phys = (p_out_phys,) + @testset "differential cross section" begin + @testset "unsafe compute" begin + PS_POINT = generate_phase_space( + TESTPROC, TESTMODEL, TESTPSDEF, p_in_phys, p_out_phys + ) - p_combs_phys = Iterators.product(p_in_all_phys, p_out_all_phys) + diffCS_on_psp = unsafe_differential_cross_section(PS_POINT) + groundtruth = TestImplementation._groundtruth_unsafe_diffCS( + TESTPROC, p_in_phys, p_out_phys + ) - @testset "cross section" begin - @testset "unsafe compute" begin - for (P_IN, P_OUT) in p_combs_phys - PS_POINT = generate_phase_space(TESTPROC, TESTMODEL, TESTPSDEF, P_IN, P_OUT) - diffCS_on_psp = unsafe_differential_cross_section(PS_POINT) - groundtruth = TestImplementation._groundtruth_unsafe_diffCS( - TESTPROC, P_IN, P_OUT - ) - @test isapprox(diffCS_on_psp, groundtruth, atol=ATOL, rtol=RTOL) - end + @test isapprox(diffCS_on_psp, groundtruth, atol=ATOL, rtol=RTOL) end @testset "safe compute" begin - for (P_IN, P_OUT) in p_combs_valid + for (P_IN, P_OUT) in p_combs PS_POINT = generate_phase_space(TESTPROC, TESTMODEL, TESTPSDEF, P_IN, P_OUT) + diffCS_on_psp = differential_cross_section(PS_POINT) groundtruth = TestImplementation._groundtruth_safe_diffCS( TESTPROC, P_IN, P_OUT ) + @test isapprox(diffCS_on_psp, groundtruth, atol=ATOL, rtol=RTOL) end end + end - @testset "total cross section" begin - @testset "compute" begin - for P_IN in (p_in_phys,) - groundtruth = TestImplementation._groundtruth_total_cross_section(P_IN) - totCS_on_moms = total_cross_section( - TESTPROC, TESTMODEL, TESTPSDEF, P_IN - ) - - COORDS_IN = TestImplementation.flat_components(P_IN) - totCS_on_coords = total_cross_section( - TESTPROC, TESTMODEL, TESTPSDEF, COORDS_IN - ) - - @test isapprox(totCS_on_moms, groundtruth, atol=ATOL, rtol=RTOL) - @test isapprox(totCS_on_coords, groundtruth, atol=ATOL, rtol=RTOL) - end - end - @testset "invalid input" begin - for P_IN in (p_in_phys_invalid,) - @test_throws DimensionMismatch total_cross_section( - TESTPROC, TESTMODEL, TESTPSDEF, P_IN - ) - - COORDS_IN = TestImplementation.flat_components(P_IN) - @test_throws DimensionMismatch total_cross_section( - TESTPROC, TESTMODEL, TESTPSDEF, COORDS_IN - ) - end - end + @testset "total cross section" begin + @testset "compute" begin + groundtruth = TestImplementation._groundtruth_total_cross_section(p_in_phys) + totCS_on_moms = total_cross_section(TESTPROC, TESTMODEL, TESTPSDEF, p_in_phys) + + COORDS_IN = TestImplementation.flat_components(p_in_phys) + totCS_on_coords = total_cross_section(TESTPROC, TESTMODEL, TESTPSDEF, COORDS_IN) + + @test isapprox(totCS_on_moms, groundtruth, atol=ATOL, rtol=RTOL) + @test isapprox(totCS_on_coords, groundtruth, atol=ATOL, rtol=RTOL) + end + @testset "invalid input" begin + @test_throws DimensionMismatch total_cross_section( + TESTPROC, TESTMODEL, TESTPSDEF, p_in_phys_invalid + ) + + COORDS_IN = TestImplementation.flat_components(p_in_phys_invalid) + @test_throws DimensionMismatch total_cross_section( + TESTPROC, TESTMODEL, TESTPSDEF, COORDS_IN + ) end + end - @testset "differential probability" begin - @testset "unsafe compute" begin - for (P_IN, P_OUT) in p_combs_phys - PS_POINT = generate_phase_space( - TESTPROC, TESTMODEL, TESTPSDEF, P_IN, P_OUT - ) - prop_on_psp = unsafe_differential_probability(PS_POINT) - groundtruth = TestImplementation._groundtruth_unsafe_probability( - TESTPROC, P_IN, P_OUT - ) - @test isapprox(prop_on_psp, groundtruth, atol=ATOL, rtol=RTOL) - end - end + @testset "differential probability" begin + @testset "unsafe compute" begin + PS_POINT = generate_phase_space( + TESTPROC, TESTMODEL, TESTPSDEF, p_in_phys, p_out_phys + ) + prop_on_psp = unsafe_differential_probability(PS_POINT) + groundtruth = TestImplementation._groundtruth_unsafe_probability( + TESTPROC, p_in_phys, p_out_phys + ) + @test isapprox(prop_on_psp, groundtruth, atol=ATOL, rtol=RTOL) + end - @testset "safe compute" begin - for (P_IN, P_OUT) in p_combs_valid - PS_POINT = generate_phase_space( - TESTPROC, TESTMODEL, TESTPSDEF, P_IN, P_OUT - ) - prop_on_psp = differential_probability(PS_POINT) - groundtruth = TestImplementation._groundtruth_safe_probability( - TESTPROC, P_IN, P_OUT - ) - @test isapprox(prop_on_psp, groundtruth, atol=ATOL, rtol=RTOL) - end + @testset "safe compute" begin + for (P_IN, P_OUT) in p_combs + PS_POINT = generate_phase_space(TESTPROC, TESTMODEL, TESTPSDEF, P_IN, P_OUT) + prop_on_psp = differential_probability(PS_POINT) + groundtruth = TestImplementation._groundtruth_safe_probability( + TESTPROC, P_IN, P_OUT + ) + @test isapprox(prop_on_psp, groundtruth, atol=ATOL, rtol=RTOL) end + end + end - @testset "total probability" begin - @testset "compute" begin - for P_IN in (p_in_phys,) - groundtruth = TestImplementation._groundtruth_total_probability( - P_IN - ) - totCS_on_moms = total_probability( - TESTPROC, TESTMODEL, TESTPSDEF, P_IN - ) - - COORDS_IN = TestImplementation.flat_components(P_IN) - totCS_on_coords = total_probability( - TESTPROC, TESTMODEL, TESTPSDEF, COORDS_IN - ) - - @test isapprox(totCS_on_moms, groundtruth, atol=ATOL, rtol=RTOL) - @test isapprox(totCS_on_coords, groundtruth, atol=ATOL, rtol=RTOL) - end - end - @testset "invalid input" begin - for P_IN in (p_in_phys_invalid,) - @test_throws DimensionMismatch total_probability( - TESTPROC, TESTMODEL, TESTPSDEF, P_IN - ) - - COORDS_IN = TestImplementation.flat_components(P_IN) - @test_throws DimensionMismatch total_probability( - TESTPROC, TESTMODEL, TESTPSDEF, COORDS_IN - ) - end - end - end + @testset "total probability" begin + @testset "compute" begin + groundtruth = TestImplementation._groundtruth_total_probability(p_in_phys) + totCS_on_moms = total_probability(TESTPROC, TESTMODEL, TESTPSDEF, p_in_phys) + + COORDS_IN = TestImplementation.flat_components(p_in_phys) + totCS_on_coords = total_probability(TESTPROC, TESTMODEL, TESTPSDEF, COORDS_IN) + + @test isapprox(totCS_on_moms, groundtruth, atol=ATOL, rtol=RTOL) + @test isapprox(totCS_on_coords, groundtruth, atol=ATOL, rtol=RTOL) + end + @testset "invalid input" begin + @test_throws DimensionMismatch total_probability( + TESTPROC, TESTMODEL, TESTPSDEF, p_in_phys_invalid + ) + + COORDS_IN = TestImplementation.flat_components(p_in_phys_invalid) + @test_throws DimensionMismatch total_probability( + TESTPROC, TESTMODEL, TESTPSDEF, COORDS_IN + ) end end end From 7e78c2b686f6d878b5ffacbb3d06010cb62a99c5 Mon Sep 17 00:00:00 2001 From: Uwe Hernandez Acosta Date: Tue, 21 May 2024 09:30:23 +0200 Subject: [PATCH 10/22] Apply suggestions from code review Co-authored-by: Anton Reinhard --- src/interfaces/process_interface.jl | 10 +++++----- .../one_photon_compton/perturbative/cross_section.jl | 7 +++---- 2 files changed, 8 insertions(+), 9 deletions(-) diff --git a/src/interfaces/process_interface.jl b/src/interfaces/process_interface.jl index 3c6c0e4..c432c1f 100644 --- a/src/interfaces/process_interface.jl +++ b/src/interfaces/process_interface.jl @@ -34,14 +34,14 @@ interface functions need to be implemented for every combination of `CustomProce _is_in_phasespace(psp::PhaseSpacePoint{CustomProcess,CustomModel}) - _phase_space_factor(psp::PhaseSpacePoint{CustomProcess,CustomModel, CustomPhasespaceDefinition}) + _phase_space_factor(psp::PhaseSpacePoint{CustomProcess,CustomModel,CustomPhasespaceDefinition}) ``` Optional is the implementation of ```Julia - _total_probability(psp::PhaseSpacePoint{CustomProcess,CustomModel, CustomPhasespaceDefinition}) + _total_probability(psp::PhaseSpacePoint{CustomProcess,CustomModel,CustomPhasespaceDefinition}) ``` to enable the calculation of total probabilities and cross sections. @@ -83,12 +83,12 @@ Interface function, which returns the incident flux of the given scattering proc function _incident_flux end """ - _matrix_elemen(PhaseSpacePoint{PROC,MODEL}) where { + _matrix_element(PhaseSpacePoint{PROC,MODEL}) where { PROC <: AbstractProcessDefinition, MODEL <: AbstractModelDefinition, } -Interface function, which returns a tuple of scattering matrix elements for each spin and polarisation combination of `proc`. +Interface function which returns a tuple of scattering matrix elements for each spin and polarisation combination of `proc`. """ function _matrix_element end @@ -108,7 +108,7 @@ function _averaging_norm end MODEL <: AbstractModelDefinition, } -Interface function, which returns `true`, if the combination of the given incoming and outgoing phase space +Interface function which returns `true` if the combination of the given incoming and outgoing phase space is physical, i.e. all momenta are on-shell and some sort of energy-momentum conservation holds. """ function _is_in_phasespace end diff --git a/src/processes/one_photon_compton/perturbative/cross_section.jl b/src/processes/one_photon_compton/perturbative/cross_section.jl index c4a5891..cb1173a 100644 --- a/src/processes/one_photon_compton/perturbative/cross_section.jl +++ b/src/processes/one_photon_compton/perturbative/cross_section.jl @@ -40,11 +40,10 @@ end function _is_in_phasespace(psp::PhaseSpacePoint{<:Compton,PerturbativeQED}) in_ps = momentum.(psp.in_particles) out_ps = momentum.(psp.out_particles) - return if (isapprox(sum(in_ps), sum(out_ps))) - _all_onshell(psp.proc, in_ps, out_ps) - else - false + if (!isapprox(sum(in_ps), sum(out_ps))) + return false end + return _all_onshell(psp.proc, in_ps, out_ps) end @inline function _phase_space_factor(psp::PhaseSpacePoint{<:Compton,PerturbativeQED}) From 94a560c3882ec1e5979f74305d648ef285116632 Mon Sep 17 00:00:00 2001 From: Anton Reinhard Date: Sun, 19 May 2024 01:36:50 +0200 Subject: [PATCH 11/22] Improve type info for PS --- src/phase_spaces.jl | 94 +++++++++++++++----------------------------- test/phase_spaces.jl | 27 +++++++------ 2 files changed, 47 insertions(+), 74 deletions(-) diff --git a/src/phase_spaces.jl b/src/phase_spaces.jl index 48e6c51..e771fd2 100644 --- a/src/phase_spaces.jl +++ b/src/phase_spaces.jl @@ -139,12 +139,9 @@ Representation of a particle with a state. It has four fields: - `dir::ParticleDirection`: The direction of the particle, `QEDbase.Incoming()` or `QEDbase.Outgoing()`. - `species::AbstractParticleType`: The species of the particle, `QEDbase.Electron()`, `QEDbase.Positron()` etc. - `mom::AbstractFourMomentum`: The momentum of the particle. -- `spin_or_pol::AbstractSpinOrPolarization`: The spin or polarization of the particle, `QEDbase.SpinUp()`, `QEDbase.PolX() etc. Can only use spins with fermions and polarizations with bosons. `AllSpin` or `AllPol` by default. Overloads for `QEDbase.is_fermion`, `QEDbase.is_boson`, `QEDbase.is_particle`, `QEDbase.is_anti_particle`, `QEDbase.is_incoming`, `QEDbase.is_outgoing`, `QEDbase.mass`, and `QEDbase.charge` are provided, delegating the call to the correct field and thus implementing the `QEDbase.AbstractParticle` interface. -The legality of the combination of `species` and `spin_or_pol` is checked on construction. If, for example, the construction of an `Electron()` with a polarization is attempted, an [`InvalidInputError`](@ref) is thrown. - ```jldoctest julia> using QEDbase; using QEDprocesses @@ -159,24 +156,19 @@ ParticleStateful: outgoing photon momentum: [1.0, 0.0, 0.0, 0.0] ``` """ -struct ParticleStateful{ElType<:AbstractFourMomentum} <: AbstractParticle - dir::ParticleDirection - species::AbstractParticleType - mom::ElType - spin_or_pol::AbstractSpinOrPolarization - - function ParticleStateful( - dir::ParticleDirection, species::Species, mom::ElType, spin::Spin=AllSpin() - ) where {Species<:FermionLike,ElType<:AbstractFourMomentum,Spin<:AbstractSpin} - # constructor for fermions with spin - return new{ElType}(dir, species, mom, spin) - end +struct ParticleStateful{ + DIR<:ParticleDirection,SPECIES<:AbstractParticleType,ELEMENT<:AbstractFourMomentum +} <: AbstractParticle + dir::DIR + species::SPECIES + mom::ELEMENT function ParticleStateful( - dir::ParticleDirection, species::Species, mom::ElType, pol::Pol=AllPol() - ) where {Species<:BosonLike,ElType<:AbstractFourMomentum,Pol<:AbstractPolarization} - # constructor for bosons with polarization - return new{ElType}(dir, species, mom, pol) + dir::DIR, species::SPECIES, mom::ELEMENT + ) where { + DIR<:ParticleDirection,SPECIES<:AbstractParticleType,ELEMENT<:AbstractFourMomentum + } + return new{DIR,SPECIES,ELEMENT}(dir, species, mom) end end @@ -195,14 +187,6 @@ particle_direction(part::ParticleStateful) = part.dir particle_species(part::ParticleStateful) = part.species momentum(part::ParticleStateful) = part.mom -@inline _spin(::Species, particle::ParticleStateful) where {Species<:FermionLike} = - particle.spin_or_pol -@inline spin(particle::ParticleStateful) = _spin(particle.species, particle) - -@inline _polarization(::Species, particle::ParticleStateful) where {Species<:BosonLike} = - particle.spin_or_pol -@inline polarization(particle::ParticleStateful) = _polarization(particle.species, particle) - function Base.show(io::IO, particle::ParticleStateful) print( io, "$(particle.dir) $(particle.species) ($(particle.spin_or_pol)): $(particle.mom)" @@ -212,10 +196,6 @@ end function Base.show(io::IO, m::MIME"text/plain", particle::ParticleStateful) println(io, "ParticleStateful: $(particle.dir) $(particle.species)") - println( - io, - " $(particle.spin_or_pol isa AbstractSpin ? "spin" : "polarization"): $(particle.spin_or_pol)", - ) println(io, " momentum: $(particle.mom)") return nothing end @@ -259,26 +239,24 @@ struct PhaseSpacePoint{ PROC<:AbstractProcessDefinition, MODEL<:AbstractModelDefinition, PSDEF<:AbstractPhasespaceDefinition, - PhaseSpaceElementType<:AbstractFourMomentum, - N_IN_PARTICLES, - N_OUT_PARTICLES, + IN_PARTICLES<:Tuple, + OUT_PARTICLES<:Tuple, } proc::PROC model::MODEL ps_def::PSDEF - in_particles::SVector{N_IN_PARTICLES,ParticleStateful{PhaseSpaceElementType}} - out_particles::SVector{N_OUT_PARTICLES,ParticleStateful{PhaseSpaceElementType}} + in_particles::IN_PARTICLES + out_particles::OUT_PARTICLES function PhaseSpacePoint( - proc::PROC, model::MODEL, ps_def::PSDEF, in_p::IN_P, out_p::OUT_P + proc::PROC, model::MODEL, ps_def::PSDEF, in_p::IN_PARTICLES, out_p::OUT_PARTICLES ) where { PROC<:AbstractProcessDefinition, MODEL<:AbstractModelDefinition, PSDEF<:AbstractPhasespaceDefinition, - PhaseSpaceElementType<:AbstractFourMomentum, - IN_P<:AbstractVector{ParticleStateful{PhaseSpaceElementType}}, - OUT_P<:AbstractVector{ParticleStateful{PhaseSpaceElementType}}, + IN_PARTICLES<:Tuple, + OUT_PARTICLES<:Tuple, } length(incoming_particles(proc)) == length(in_p) || throw( InvalidInputError( @@ -316,7 +294,7 @@ struct PhaseSpacePoint{ ) end - return new{PROC,MODEL,PSDEF,PhaseSpaceElementType,length(in_p),length(out_p)}( + return new{PROC,MODEL,PSDEF,IN_PARTICLES,OUT_PARTICLES}( proc, model, ps_def, in_p, out_p ) end @@ -354,9 +332,9 @@ end proc::AbstractProcessDefinition, model::AbstractModelDefinition, ps_def::AbstractPhasespaceDefinition, - in_ps::AbstractVector{ElType}, - out_ps::AbstractVector{ElType}, - ) where {ElType<:QEDbase.AbstractFourMomentum} + in_ps::AbstractVector{ELEMENT}, + out_ps::AbstractVector{ELEMENT}, + ) where {ELEMENT<:QEDbase.AbstractFourMomentum} Return the respective phase space point for given momenta of incoming and outgoing particles regarding a given process. """ @@ -364,28 +342,20 @@ function generate_phase_space( proc::AbstractProcessDefinition, model::AbstractModelDefinition, ps_def::AbstractPhasespaceDefinition, - in_ps::AbstractVector{ElType}, - out_ps::AbstractVector{ElType}, -) where {ElType<:QEDbase.AbstractFourMomentum} - in_particles = incoming_particles(proc) - in_n = number_incoming_particles(proc) - in_parts = SVector{in_n,ParticleStateful{SFourMomentum}}( - collect( - ParticleStateful(Incoming(), particle, mom) for - (particle, mom) in zip(in_particles, in_ps) - ), + in_ps::AbstractVector{ELEMENT}, + out_ps::AbstractVector{ELEMENT}, +) where {ELEMENT<:QEDbase.AbstractFourMomentum} + in_particles = Tuple( + ParticleStateful(Incoming(), particle, mom) for + (particle, mom) in zip(incoming_particles(proc), in_ps) ) - out_particles = outgoing_particles(proc) - out_n = number_outgoing_particles(proc) - out_parts = SVector{out_n,ParticleStateful{SFourMomentum}}( - collect( - ParticleStateful(Outgoing(), particle, mom) for - (particle, mom) in zip(out_particles, out_ps) - ), + out_particles = Tuple( + ParticleStateful(Outgoing(), particle, mom) for + (particle, mom) in zip(outgoing_particles(proc), out_ps) ) - return PhaseSpacePoint(proc, model, ps_def, in_parts, out_parts) + return PhaseSpacePoint(proc, model, ps_def, in_particles, out_particles) end function Base.show(io::IO, psp::PhaseSpacePoint) diff --git a/test/phase_spaces.jl b/test/phase_spaces.jl index 817d6d3..c4f97fc 100644 --- a/test/phase_spaces.jl +++ b/test/phase_spaces.jl @@ -21,7 +21,6 @@ end @testset "Stateful Particle" begin DIRECTIONS = [Incoming(), Outgoing()] SPECIES = [Electron(), Positron()] #=, Muon(), AntiMuon(), Tauon(), AntiTauon()=# - SPINANDPOLS = [AllSpin(), SpinUp(), SpinDown(), AllPol(), PolX(), PolY()] for (species, dir, spin_or_pol) in Iterators.product(SPECIES, DIRECTIONS, SPINANDPOLS) mom = rand(RNG, SFourMomentum) @@ -85,11 +84,11 @@ end out_el = ParticleStateful(Outgoing(), Electron(), out_el_mom) out_ph = ParticleStateful(Outgoing(), Photon(), out_ph_mom) - in_particles_valid = SVector(in_el, in_ph) - in_particles_invalid = SVector(in_el, out_ph) + in_particles_valid = (in_el, in_ph) + in_particles_invalid = (in_el, out_ph) - out_particles_valid = SVector(out_el, out_ph) - out_particles_invalid = SVector(out_el, in_ph) + out_particles_valid = (out_el, out_ph) + out_particles_invalid = (out_el, in_ph) model = TESTMODEL process = TestImplementation.TestProcess( @@ -134,12 +133,12 @@ end process, model, phasespace_def, in_particles_valid, out_particles_invalid ) - @test_throws r"process given particle species \(electron\) does not match stateful particle species \(photon\)" PhaseSpacePoint( - process, model, phasespace_def, SVector(in_ph, in_el), out_particles_valid + @test_throws r"process given particle species \((.*)Electron\(\)\) does not match stateful particle species \((.*)Photon\(\)\)" PhaseSpacePoint( + process, model, phasespace_def, (in_ph, in_el), out_particles_valid ) - @test_throws r"process given particle species \(electron\) does not match stateful particle species \(photon\)" PhaseSpacePoint( - process, model, phasespace_def, in_particles_valid, SVector(out_ph, out_el) + @test_throws r"process given particle species \((.*)Electron\(\)\) does not match stateful particle species \((.*)Photon\(\)\)" PhaseSpacePoint( + process, model, phasespace_def, in_particles_valid, (out_ph, out_el) ) end @@ -162,17 +161,21 @@ end ) @test_throws InvalidInputError PhaseSpacePoint( - process, model, phasespace_def, SVector(in_ph, in_el), out_particles_valid + process, model, phasespace_def, (in_ph, in_el), out_particles_valid ) @test_throws InvalidInputError PhaseSpacePoint( - process, model, phasespace_def, in_particles_valid, SVector(out_ph, out_el) + process, model, phasespace_def, in_particles_valid, (out_ph, out_el) ) end @testset "Generation" begin test_psp = generate_phase_space( - process, model, phasespace_def, [in_el_mom, in_ph_mom], [out_el_mom, out_ph_mom] + process, + model, + phasespace_def, + SVector(in_el_mom, in_ph_mom), + SVector(out_el_mom, out_ph_mom), ) @test test_psp.proc == process From f2b646c3ae7d615e2fafed1f769c811e269fec29 Mon Sep 17 00:00:00 2001 From: Anton Reinhard Date: Sun, 19 May 2024 01:39:53 +0200 Subject: [PATCH 12/22] Move generate_phasespacepoint into constructor --- src/phase_spaces.jl | 62 ++++++++++++++++++++++---------------------- test/phase_spaces.jl | 2 +- 2 files changed, 32 insertions(+), 32 deletions(-) diff --git a/src/phase_spaces.jl b/src/phase_spaces.jl index e771fd2..d0071b8 100644 --- a/src/phase_spaces.jl +++ b/src/phase_spaces.jl @@ -298,6 +298,37 @@ struct PhaseSpacePoint{ proc, model, ps_def, in_p, out_p ) end + + """ + PhaseSpacePoint( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + ps_def::AbstractPhasespaceDefinition, + in_ps::AbstractVector{ELEMENT}, + out_ps::AbstractVector{ELEMENT}, + ) where {ELEMENT<:QEDbase.AbstractFourMomentum} + + Construct the phase space point from given momenta of incoming and outgoing particles regarding a given process. + """ + function PhaseSpacePoint( + proc::AbstractProcessDefinition, + model::AbstractModelDefinition, + ps_def::AbstractPhasespaceDefinition, + in_ps::AbstractVector{ELEMENT}, + out_ps::AbstractVector{ELEMENT}, + ) where {ELEMENT<:QEDbase.AbstractFourMomentum} + in_particles = Tuple( + ParticleStateful(Incoming(), particle, mom) for + (particle, mom) in zip(incoming_particles(proc), in_ps) + ) + + out_particles = Tuple( + ParticleStateful(Outgoing(), particle, mom) for + (particle, mom) in zip(outgoing_particles(proc), out_ps) + ) + + return PhaseSpacePoint(proc, model, ps_def, in_particles, out_particles) + end end """ @@ -327,37 +358,6 @@ function momentum(psp::PhaseSpacePoint, dir::ParticleDirection, n::Int) return psp[dir, n].mom end -""" - generate_phase_space( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - ps_def::AbstractPhasespaceDefinition, - in_ps::AbstractVector{ELEMENT}, - out_ps::AbstractVector{ELEMENT}, - ) where {ELEMENT<:QEDbase.AbstractFourMomentum} - -Return the respective phase space point for given momenta of incoming and outgoing particles regarding a given process. -""" -function generate_phase_space( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - ps_def::AbstractPhasespaceDefinition, - in_ps::AbstractVector{ELEMENT}, - out_ps::AbstractVector{ELEMENT}, -) where {ELEMENT<:QEDbase.AbstractFourMomentum} - in_particles = Tuple( - ParticleStateful(Incoming(), particle, mom) for - (particle, mom) in zip(incoming_particles(proc), in_ps) - ) - - out_particles = Tuple( - ParticleStateful(Outgoing(), particle, mom) for - (particle, mom) in zip(outgoing_particles(proc), out_ps) - ) - - return PhaseSpacePoint(proc, model, ps_def, in_particles, out_particles) -end - function Base.show(io::IO, psp::PhaseSpacePoint) print(io, "PhaseSpacePoint of $(psp.proc)") return nothing diff --git a/test/phase_spaces.jl b/test/phase_spaces.jl index c4f97fc..25a3c28 100644 --- a/test/phase_spaces.jl +++ b/test/phase_spaces.jl @@ -170,7 +170,7 @@ end end @testset "Generation" begin - test_psp = generate_phase_space( + test_psp = PhaseSpacePoint( process, model, phasespace_def, From a00d35ba6e1a957333e7725370a15308733acdfb Mon Sep 17 00:00:00 2001 From: Anton Reinhard Date: Sun, 19 May 2024 05:40:54 +0200 Subject: [PATCH 13/22] Improve performance with type magic --- src/interfaces/process_interface.jl | 9 ++- src/phase_spaces.jl | 111 ++++++++++++++++------------ 2 files changed, 72 insertions(+), 48 deletions(-) diff --git a/src/interfaces/process_interface.jl b/src/interfaces/process_interface.jl index 8400ac4..c889306 100644 --- a/src/interfaces/process_interface.jl +++ b/src/interfaces/process_interface.jl @@ -180,7 +180,6 @@ function _phase_space_factor end ####################### """ - number_incoming_particles(proc_def::AbstractProcessDefinition) Return the number of incoming particles of a given process. @@ -190,7 +189,6 @@ Return the number of incoming particles of a given process. end """ - number_outgoing_particles(proc_def::AbstractProcessDefinition) Return the number of outgoing particles of a given process. @@ -199,6 +197,13 @@ Return the number of outgoing particles of a given process. return length(outgoing_particles(proc_def)) end +@inline function incoming_particle_types(proc_def::AbstractProcessDefinition)::Type + return typeof(incoming_particles(proc_def)) +end +@inline function outgoing_particle_types(proc_def::AbstractProcessDefinition)::Type + return typeof(outgoing_particles(proc_def)) +end + """ in_phase_space_dimension( proc::AbstractProcessDefinition, diff --git a/src/phase_spaces.jl b/src/phase_spaces.jl index d0071b8..eacebc4 100644 --- a/src/phase_spaces.jl +++ b/src/phase_spaces.jl @@ -239,8 +239,8 @@ struct PhaseSpacePoint{ PROC<:AbstractProcessDefinition, MODEL<:AbstractModelDefinition, PSDEF<:AbstractPhasespaceDefinition, - IN_PARTICLES<:Tuple, - OUT_PARTICLES<:Tuple, + IN_PARTICLES<:Tuple{Vararg{ParticleStateful}}, + OUT_PARTICLES<:Tuple{Vararg{ParticleStateful}}, } proc::PROC model::MODEL @@ -255,44 +255,10 @@ struct PhaseSpacePoint{ PROC<:AbstractProcessDefinition, MODEL<:AbstractModelDefinition, PSDEF<:AbstractPhasespaceDefinition, - IN_PARTICLES<:Tuple, - OUT_PARTICLES<:Tuple, + IN_PARTICLES<:Tuple{Vararg{ParticleStateful}}, + OUT_PARTICLES<:Tuple{Vararg{ParticleStateful}}, } - length(incoming_particles(proc)) == length(in_p) || throw( - InvalidInputError( - "the number of incoming particles given by the process ($(incoming_particles(proc))) mismatches the number of given stateful incoming particles ($(length(in_p)))", - ), - ) - length(outgoing_particles(proc)) == length(out_p) || throw( - InvalidInputError( - "the number of outgoing particles given by the process ($(outgoing_particles(proc))) mismatches the number of given stateful outgoing particles ($(length(out_p)))", - ), - ) - - for (proc_p, p) in zip(incoming_particles(proc), in_p) - proc_p == p.species || throw( - InvalidInputError( - "process given particle species ($(proc_p)) does not match stateful particle species ($(p.species))", - ), - ) - is_incoming(p) || throw( - InvalidInputError( - "stateful particle $(p) is given as an incoming particle but is outgoing", - ), - ) - end - for (proc_p, p) in zip(outgoing_particles(proc), out_p) - proc_p == p.species || throw( - InvalidInputError( - "process given particle species ($(proc_p)) does not match stateful particle species ($(p.species))", - ), - ) - is_outgoing(p) || throw( - InvalidInputError( - "stateful particle $(p) is given as an outgoing particle but is incoming", - ), - ) - end + _check_psp(incoming_particles(proc), outgoing_particles(proc), in_p, out_p) return new{PROC,MODEL,PSDEF,IN_PARTICLES,OUT_PARTICLES}( proc, model, ps_def, in_p, out_p @@ -311,26 +277,79 @@ struct PhaseSpacePoint{ Construct the phase space point from given momenta of incoming and outgoing particles regarding a given process. """ function PhaseSpacePoint( - proc::AbstractProcessDefinition, - model::AbstractModelDefinition, - ps_def::AbstractPhasespaceDefinition, + proc::PROC, + model::MODEL, + ps_def::PSDEF, in_ps::AbstractVector{ELEMENT}, out_ps::AbstractVector{ELEMENT}, - ) where {ELEMENT<:QEDbase.AbstractFourMomentum} - in_particles = Tuple( + ) where { + PROC<:AbstractProcessDefinition, + MODEL<:AbstractModelDefinition, + PSDEF<:AbstractPhasespaceDefinition, + ELEMENT<:AbstractFourMomentum, + } + in_particles = Tuple{ + Vararg{ParticleStateful{Incoming},number_incoming_particles(proc)} + }( ParticleStateful(Incoming(), particle, mom) for (particle, mom) in zip(incoming_particles(proc), in_ps) ) - out_particles = Tuple( + out_particles = Tuple{ + Vararg{ParticleStateful{Outgoing},number_outgoing_particles(proc)} + }( ParticleStateful(Outgoing(), particle, mom) for (particle, mom) in zip(outgoing_particles(proc), out_ps) ) - return PhaseSpacePoint(proc, model, ps_def, in_particles, out_particles) + # no need to check anything since it is constructed correctly above + + return new{PROC,MODEL,PSDEF,typeof(in_particles),typeof(out_particles)}( + proc, model, ps_def, in_particles, out_particles + ) end end +@inline function _single_type_check(t::ParticleStateful, p::AbstractParticleType) + t.species == p || throw( + InvalidInputError( + "process given particle species ($(p)) does not match stateful particle species ($(t.species))", + ), + ) + return nothing +end + +@inline function _check_psp( + in_proc::P_IN_Ts, out_proc::P_OUT_Ts, in_p::IN_Ts, out_p::OUT_Ts +) where { + P_IN_Ts<:Tuple{Vararg{AbstractParticleType}}, + P_OUT_Ts<:Tuple{Vararg{AbstractParticleType}}, + IN_Ts<:Tuple{Vararg{ParticleStateful{Incoming}}}, + OUT_Ts<:Tuple{Vararg{ParticleStateful{Outgoing}}}, +} + length(in_proc) == length(in_p) || throw( + InvalidInputError( + "the number of incoming particles given by the process ($(length(in_proc))) mismatches the number of given stateful incoming particles ($(length(in_p)))", + ), + ) + length(out_proc) == length(out_p) || throw( + InvalidInputError( + "the number of outgoing particles given by the process ($(length(out_proc))) mismatches the number of given stateful outgoing particles ($(length(out_p)))", + ), + ) + + # in_proc/out_proc contain only species types + # in_p/out_p contain full ParticleStateful types + + for (in_species_p, t) in zip(in_proc, in_p) + _single_type_check(t, in_species_p) + end + for (out_species_p, t) in zip(out_proc, out_p) + _single_type_check(t, out_species_p) + end + return nothing +end + """ Base.getindex(psp::PhaseSpacePoint, dir::Incoming, n::Int) From 7a1fa95e888795c761b9b324c26c9e9ff78d7c09 Mon Sep 17 00:00:00 2001 From: AntonReinhard Date: Sun, 19 May 2024 05:50:22 +0200 Subject: [PATCH 14/22] EOD --- test/cross_sections.jl | 8 ++++---- test/phase_spaces.jl | 31 +++++++++++-------------------- 2 files changed, 15 insertions(+), 24 deletions(-) diff --git a/test/cross_sections.jl b/test/cross_sections.jl index 7c26a6f..f9127c2 100644 --- a/test/cross_sections.jl +++ b/test/cross_sections.jl @@ -119,7 +119,7 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() end @testset "compute on phase space points" begin - PS_POINT = generate_phase_space( + PS_POINT = PhaseSpacePoint( TESTPROC, TESTMODEL, TESTPSDEF, p_in_phys, p_out_phys ) diffCS_on_psp = unsafe_differential_cross_section(PS_POINT) @@ -168,7 +168,7 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() end @testset "compute on phase space points" begin - PS_POINT = generate_phase_space( + PS_POINT = PhaseSpacePoint( TESTPROC, TESTMODEL, TESTPSDEF, p_in_phys, p_out_phys ) diffCS_on_psp = differential_cross_section(PS_POINT) @@ -252,7 +252,7 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() end @testset "compute on phase space points" begin - PS_POINT = generate_phase_space( + PS_POINT = PhaseSpacePoint( TESTPROC, TESTMODEL, TESTPSDEF, p_in_phys, p_out_phys ) prop_on_psp = unsafe_differential_probability(PS_POINT) @@ -301,7 +301,7 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() end @testset "compute on phase space points" begin - PS_POINT = generate_phase_space( + PS_POINT = PhaseSpacePoint( TESTPROC, TESTMODEL, TESTPSDEF, p_in_phys, p_out_phys ) prop_on_psp = differential_probability(PS_POINT) diff --git a/test/phase_spaces.jl b/test/phase_spaces.jl index 25a3c28..2c427b8 100644 --- a/test/phase_spaces.jl +++ b/test/phase_spaces.jl @@ -22,22 +22,20 @@ end DIRECTIONS = [Incoming(), Outgoing()] SPECIES = [Electron(), Positron()] #=, Muon(), AntiMuon(), Tauon(), AntiTauon()=# - for (species, dir, spin_or_pol) in Iterators.product(SPECIES, DIRECTIONS, SPINANDPOLS) + for (species, dir) in Iterators.product(SPECIES, DIRECTIONS) mom = rand(RNG, SFourMomentum) - if (is_fermion(species) && (spin_or_pol isa AbstractSpin)) || - (is_boson(species) && (spin_or_pol isa AbstractPolarization)) - particle_stateful = ParticleStateful(dir, species, mom, spin_or_pol) + particle_stateful = ParticleStateful(dir, species, mom) - # particle interface - @test is_fermion(particle_stateful) == is_fermion(species) - @test is_boson(particle_stateful) == is_boson(species) - @test is_particle(particle_stateful) == is_particle(species) - @test is_anti_particle(particle_stateful) == is_anti_particle(species) - @test is_incoming(particle_stateful) == is_incoming(dir) - @test is_outgoing(particle_stateful) == is_outgoing(dir) - @test mass(particle_stateful) == mass(species) - @test charge(particle_stateful) == charge(species) + # particle interface + @test is_fermion(particle_stateful) == is_fermion(species) + @test is_boson(particle_stateful) == is_boson(species) + @test is_particle(particle_stateful) == is_particle(species) + @test is_anti_particle(particle_stateful) == is_anti_particle(species) + @test is_incoming(particle_stateful) == is_incoming(dir) + @test is_outgoing(particle_stateful) == is_outgoing(dir) + @test mass(particle_stateful) == mass(species) + @test charge(particle_stateful) == charge(species) # accessors @test particle_stateful.dir == dir @@ -46,13 +44,6 @@ end @test particle_species(particle_stateful) == particle_stateful.species @test particle_stateful.mom == mom @test momentum(particle_stateful) == mom - if (is_fermion(species)) - @test spin(particle_stateful) == spin_or_pol - @test_throws MethodError polarization(particle_stateful) - else - @test polarization(particle_stateful) == spin_or_pol - @test_throws MethodError spin(particle_stateful) - end # printing print(BUF, particle_stateful) From 1010dae3d136f77d8d2bbd1a88147692a87f3a6c Mon Sep 17 00:00:00 2001 From: Anton Reinhard Date: Mon, 20 May 2024 12:18:20 +0200 Subject: [PATCH 15/22] Improve type checking with even more type magic and fix (some) tests --- src/phase_spaces.jl | 66 +++++++++++++++--------- test/cross_sections.jl | 4 +- test/interfaces/process_interface.jl | 4 +- test/phase_spaces.jl | 13 ++--- test/test_implementation/test_process.jl | 8 +-- 5 files changed, 56 insertions(+), 39 deletions(-) diff --git a/src/phase_spaces.jl b/src/phase_spaces.jl index eacebc4..62c32fa 100644 --- a/src/phase_spaces.jl +++ b/src/phase_spaces.jl @@ -310,43 +310,63 @@ struct PhaseSpacePoint{ end end -@inline function _single_type_check(t::ParticleStateful, p::AbstractParticleType) - t.species == p || throw( +# recursion termination: success +@inline _recursive_type_check(t::Tuple{}, p::Tuple{}, dir::ParticleDirection) = nothing + +# recursion termination: overload for unequal number of particles +@inline function _recursive_type_check( + ::Tuple{Vararg{ParticleStateful,N}}, + ::Tuple{Vararg{AbstractParticleType,M}}, + dir::ParticleDirection, +) where {N,M} + return throw( InvalidInputError( - "process given particle species ($(p)) does not match stateful particle species ($(t.species))", + "the number of $(dir) particles in the process $(M) does not match number of particles in the input $(N)", ), ) - return nothing end -@inline function _check_psp( - in_proc::P_IN_Ts, out_proc::P_OUT_Ts, in_p::IN_Ts, out_p::OUT_Ts +# recursion termination: overload for invalid types +@inline function _recursive_type_check( + ::Tuple{ParticleStateful{DIR_IN_T,SPECIES_IN_T},Vararg{ParticleStateful,N}}, + ::Tuple{SPECIES_T,Vararg{AbstractParticleType,M}}, + dir::DIR_T, ) where { - P_IN_Ts<:Tuple{Vararg{AbstractParticleType}}, - P_OUT_Ts<:Tuple{Vararg{AbstractParticleType}}, - IN_Ts<:Tuple{Vararg{ParticleStateful{Incoming}}}, - OUT_Ts<:Tuple{Vararg{ParticleStateful{Outgoing}}}, + N, + M, + DIR_IN_T<:ParticleDirection, + DIR_T<:ParticleDirection, + SPECIES_IN_T<:AbstractParticleType, + SPECIES_T<:AbstractParticleType, } - length(in_proc) == length(in_p) || throw( + return throw( InvalidInputError( - "the number of incoming particles given by the process ($(length(in_proc))) mismatches the number of given stateful incoming particles ($(length(in_p)))", - ), - ) - length(out_proc) == length(out_p) || throw( - InvalidInputError( - "the number of outgoing particles given by the process ($(length(out_proc))) mismatches the number of given stateful outgoing particles ($(length(out_p)))", + "expected $(dir) $(SPECIES_T()) but got $(DIR_IN_T()) $(SPECIES_IN_T())" ), ) +end +@inline function _recursive_type_check( + t::Tuple{ParticleStateful{DIR_T,SPECIES_T},Vararg{ParticleStateful,N}}, + p::Tuple{SPECIES_T,Vararg{AbstractParticleType,N}}, + dir::DIR_T, +) where {N,DIR_T<:ParticleDirection,SPECIES_T<:AbstractParticleType} + return _recursive_type_check(t[2:end], p[2:end], dir) +end + +@inline function _check_psp( + in_proc::P_IN_Ts, out_proc::P_OUT_Ts, in_p::IN_Ts, out_p::OUT_Ts +) where { + P_IN_Ts<:Tuple{Vararg{AbstractParticleType}}, + P_OUT_Ts<:Tuple{Vararg{AbstractParticleType}}, + IN_Ts<:Tuple{Vararg{ParticleStateful}}, + OUT_Ts<:Tuple{Vararg{ParticleStateful}}, +} # in_proc/out_proc contain only species types # in_p/out_p contain full ParticleStateful types - for (in_species_p, t) in zip(in_proc, in_p) - _single_type_check(t, in_species_p) - end - for (out_species_p, t) in zip(out_proc, out_p) - _single_type_check(t, out_species_p) - end + _recursive_type_check(in_p, in_proc, Incoming()) + _recursive_type_check(out_p, out_proc, Outgoing()) return nothing end diff --git a/test/cross_sections.jl b/test/cross_sections.jl index f9127c2..5cff2bb 100644 --- a/test/cross_sections.jl +++ b/test/cross_sections.jl @@ -13,8 +13,8 @@ TESTPSDEF = TestImplementation.TestPhasespaceDef() @testset "($N_INCOMING,$N_OUTGOING)" for (N_INCOMING, N_OUTGOING) in Iterators.product( (1, rand(RNG, 2:8)), (1, rand(RNG, 2:8)) ) - INCOMING_PARTICLES = rand(RNG, TestImplementation.PARTICLE_SET, N_INCOMING) - OUTGOING_PARTICLES = rand(RNG, TestImplementation.PARTICLE_SET, N_OUTGOING) + INCOMING_PARTICLES = Tuple(rand(RNG, TestImplementation.PARTICLE_SET, N_INCOMING)) + OUTGOING_PARTICLES = Tuple(rand(RNG, TestImplementation.PARTICLE_SET, N_OUTGOING)) TESTPROC = TestImplementation.TestProcess(INCOMING_PARTICLES, OUTGOING_PARTICLES) diff --git a/test/interfaces/process_interface.jl b/test/interfaces/process_interface.jl index 1936195..8ffb684 100644 --- a/test/interfaces/process_interface.jl +++ b/test/interfaces/process_interface.jl @@ -11,8 +11,8 @@ include("../test_implementation/TestImplementation.jl") @testset "($N_INCOMING,$N_OUTGOING)" for (N_INCOMING, N_OUTGOING) in Iterators.product( (1, rand(RNG, 2:8)), (1, rand(RNG, 2:8)) ) - INCOMING_PARTICLES = rand(RNG, TestImplementation.PARTICLE_SET, N_INCOMING) - OUTGOING_PARTICLES = rand(RNG, TestImplementation.PARTICLE_SET, N_OUTGOING) + INCOMING_PARTICLES = Tuple(rand(RNG, TestImplementation.PARTICLE_SET, N_INCOMING)) + OUTGOING_PARTICLES = Tuple(rand(RNG, TestImplementation.PARTICLE_SET, N_OUTGOING)) TESTPROC = TestImplementation.TestProcess(INCOMING_PARTICLES, OUTGOING_PARTICLES) TESTPROC_FAIL = TestImplementation.TestProcess_FAIL( diff --git a/test/phase_spaces.jl b/test/phase_spaces.jl index 2c427b8..89db929 100644 --- a/test/phase_spaces.jl +++ b/test/phase_spaces.jl @@ -82,10 +82,7 @@ end out_particles_invalid = (out_el, in_ph) model = TESTMODEL - process = TestImplementation.TestProcess( - SVector{2,AbstractParticle}(Electron(), Photon()), - SVector{2,AbstractParticle}(Electron(), Photon()), - ) + process = TestImplementation.TestProcess((Electron(), Photon()), (Electron(), Photon())) phasespace_def = TESTPSDEF psp = PhaseSpacePoint( @@ -116,19 +113,19 @@ end @testset "Error handling" begin if (VERSION >= v"1.8") # julia versions before 1.8 did not have support for regex matching in @test_throws - @test_throws r"stateful particle (.*) is given as an incoming particle but is outgoing" PhaseSpacePoint( + @test_throws r"expected (.*)[Ii]ncoming(.*) (.*)[Pp]hoton(.*) but got (.*)[Oo]utgoing(.*) (.*)[Pp]hoton(.*)" PhaseSpacePoint( process, model, phasespace_def, in_particles_invalid, out_particles_valid ) - @test_throws r"stateful particle (.*) is given as an outgoing particle but is incoming" PhaseSpacePoint( + @test_throws r"expected (.*)[Oo]utgoing(.*) (.*)[Pp]hoton(.*) but got (.*)[Ii]ncoming(.*) (.*)[Pp]hoton(.*)" PhaseSpacePoint( process, model, phasespace_def, in_particles_valid, out_particles_invalid ) - @test_throws r"process given particle species \((.*)Electron\(\)\) does not match stateful particle species \((.*)Photon\(\)\)" PhaseSpacePoint( + @test_throws r"expected (.*)[Ii]ncoming(.*) (.*)[Ee]lectron(.*) but got (.*)[Ii]ncoming(.*) (.*)[Pp]hoton(.*)" PhaseSpacePoint( process, model, phasespace_def, (in_ph, in_el), out_particles_valid ) - @test_throws r"process given particle species \((.*)Electron\(\)\) does not match stateful particle species \((.*)Photon\(\)\)" PhaseSpacePoint( + @test_throws r"expected (.*)[Oo]utgoing(.*) (.*)[Ee]lectron(.*) but got (.*)[Oo]utgoing(.*) (.*)[Pp]hoton(.*)" PhaseSpacePoint( process, model, phasespace_def, in_particles_valid, (out_ph, out_el) ) end diff --git a/test/test_implementation/test_process.jl b/test/test_implementation/test_process.jl index fb7722d..58829ad 100644 --- a/test/test_implementation/test_process.jl +++ b/test/test_implementation/test_process.jl @@ -10,7 +10,7 @@ const PARTICLE_SET = [TestParticleFermion(), TestParticleBoson()] TestProcess(rng,incoming_particles,outgoing_particles) """ -struct TestProcess{IP<:AbstractVector,OP<:AbstractVector} <: AbstractProcessDefinition +struct TestProcess{IP<:Tuple,OP<:Tuple} <: AbstractProcessDefinition incoming_particles::IP outgoing_particles::OP end @@ -24,14 +24,14 @@ end QEDprocesses.incoming_particles(proc::TestProcess) = proc.incoming_particles QEDprocesses.outgoing_particles(proc::TestProcess) = proc.outgoing_particles -struct TestProcess_FAIL{IP<:AbstractVector,OP<:AbstractVector} <: AbstractProcessDefinition +struct TestProcess_FAIL{IP<:Tuple,OP<:Tuple} <: AbstractProcessDefinition incoming_particles::IP outgoing_particles::OP end function TestProcess_FAIL(rng::AbstractRNG, N_in::Int, N_out::Int) - in_particles = rand(rng, PARTICLE_SET, N_in) - out_particles = rand(rng, PARTICLE_SET, N_out) + in_particles = Tuple(rand(rng, PARTICLE_SET, N_in)) + out_particles = Tuple(rand(rng, PARTICLE_SET, N_out)) return TestProcess_FAIL(in_particles, out_particles) end From 2cf05e607e795994899e2782e1675a318db20562 Mon Sep 17 00:00:00 2001 From: Anton Reinhard Date: Mon, 20 May 2024 12:36:43 +0200 Subject: [PATCH 16/22] Fix last type instability in psp construction from momenta --- src/phase_spaces.jl | 17 +++++++++++++++-- 1 file changed, 15 insertions(+), 2 deletions(-) diff --git a/src/phase_spaces.jl b/src/phase_spaces.jl index 62c32fa..34f8233 100644 --- a/src/phase_spaces.jl +++ b/src/phase_spaces.jl @@ -187,6 +187,19 @@ particle_direction(part::ParticleStateful) = part.dir particle_species(part::ParticleStateful) = part.species momentum(part::ParticleStateful) = part.mom +# recursion termination: base case +@inline _assemble_tuple_type(proc::Tuple{}, dir::ParticleDirection, ELTYPE::Type) = () + +# function assembling the correct type information for the tuple of ParticleStatefuls in a phasespace point constructed from momenta +@inline function _assemble_tuple_type( + proc::Tuple{SPECIES_T,Vararg{AbstractParticleType}}, dir::DIR_T, ELTYPE::Type +) where {SPECIES_T<:AbstractParticleType,DIR_T<:ParticleDirection} + return ( + ParticleStateful{DIR_T,SPECIES_T,ELTYPE}, + _assemble_tuple_type(proc[2:end], dir, ELTYPE)..., + ) +end + function Base.show(io::IO, particle::ParticleStateful) print( io, "$(particle.dir) $(particle.species) ($(particle.spin_or_pol)): $(particle.mom)" @@ -289,14 +302,14 @@ struct PhaseSpacePoint{ ELEMENT<:AbstractFourMomentum, } in_particles = Tuple{ - Vararg{ParticleStateful{Incoming},number_incoming_particles(proc)} + _assemble_tuple_type(incoming_particles(proc), Incoming(), ELEMENT)... }( ParticleStateful(Incoming(), particle, mom) for (particle, mom) in zip(incoming_particles(proc), in_ps) ) out_particles = Tuple{ - Vararg{ParticleStateful{Outgoing},number_outgoing_particles(proc)} + _assemble_tuple_type(outgoing_particles(proc), Outgoing(), ELEMENT)... }( ParticleStateful(Outgoing(), particle, mom) for (particle, mom) in zip(outgoing_particles(proc), out_ps) From 11250f31d747835813f13f0c5a36b077fea29fde Mon Sep 17 00:00:00 2001 From: Anton Reinhard Date: Tue, 21 May 2024 10:04:07 +0200 Subject: [PATCH 17/22] Rebase fix --- src/phase_spaces.jl | 4 +--- test/phase_spaces.jl | 42 +++++++++++++++++------------------------- 2 files changed, 18 insertions(+), 28 deletions(-) diff --git a/src/phase_spaces.jl b/src/phase_spaces.jl index 34f8233..10590e2 100644 --- a/src/phase_spaces.jl +++ b/src/phase_spaces.jl @@ -201,9 +201,7 @@ momentum(part::ParticleStateful) = part.mom end function Base.show(io::IO, particle::ParticleStateful) - print( - io, "$(particle.dir) $(particle.species) ($(particle.spin_or_pol)): $(particle.mom)" - ) + print(io, "$(particle.dir) $(particle.species): $(particle.mom)") return nothing end diff --git a/test/phase_spaces.jl b/test/phase_spaces.jl index 89db929..cd6d04f 100644 --- a/test/phase_spaces.jl +++ b/test/phase_spaces.jl @@ -37,30 +37,21 @@ end @test mass(particle_stateful) == mass(species) @test charge(particle_stateful) == charge(species) - # accessors - @test particle_stateful.dir == dir - @test particle_direction(particle_stateful) == particle_stateful.dir - @test particle_stateful.species == species - @test particle_species(particle_stateful) == particle_stateful.species - @test particle_stateful.mom == mom - @test momentum(particle_stateful) == mom - - # printing - print(BUF, particle_stateful) - @test String(take!(BUF)) == "$(dir) $(species) ($(spin_or_pol)): $(mom)" - - show(BUF, MIME"text/plain"(), particle_stateful) - @test String(take!(BUF)) == - "ParticleStateful: $(dir) $(species)\n $(spin_or_pol isa AbstractSpin ? "spin" : "polarization"): $(spin_or_pol)\n momentum: $(mom)\n" - else - if (VERSION >= v"1.8") - # julia versions before 1.8 did not have support for regex matching in @test_throws - @test_throws "MethodError: no method matching ParticleStateful" ParticleStateful( - dir, species, mom, spin_or_pol - ) - end - @test_throws MethodError ParticleStateful(dir, species, mom, spin_or_pol) - end + # accessors + @test particle_stateful.dir == dir + @test particle_direction(particle_stateful) == particle_stateful.dir + @test particle_stateful.species == species + @test particle_species(particle_stateful) == particle_stateful.species + @test particle_stateful.mom == mom + @test momentum(particle_stateful) == mom + + # printing + print(BUF, particle_stateful) + @test String(take!(BUF)) == "$(dir) $(species): $(mom)" + + show(BUF, MIME"text/plain"(), particle_stateful) + @test String(take!(BUF)) == + "ParticleStateful: $(dir) $(species)\n momentum: $(mom)\n" end end @@ -89,12 +80,13 @@ end process, model, phasespace_def, in_particles_valid, out_particles_valid ) + take!(BUF) print(BUF, psp) @test String(take!(BUF)) == "PhaseSpacePoint of $(process)" show(BUF, MIME"text/plain"(), psp) @test match( - r"PhaseSpacePoint:\n process: (.*)TestProcess(.*)\n model: (.*)TestModel(.*)\n phasespace definition: (.*)TestPhasespaceDef(.*)\n incoming particles:\n -> incoming electron \(all spins\): (.*)\n -> incoming photon \(all polarizations\): (.*)\n outgoing particles:\n -> outgoing electron \(all spins\): (.*)\n -> outgoing photon \(all polarizations\): (.*)\n", + r"PhaseSpacePoint:\n process: (.*)TestProcess(.*)\n model: (.*)TestModel(.*)\n phasespace definition: (.*)TestPhasespaceDef(.*)\n incoming particles:\n -> incoming electron: (.*)\n -> incoming photon: (.*)\n outgoing particles:\n -> outgoing electron: (.*)\n -> outgoing photon: (.*)\n", String(take!(BUF)), ) isa RegexMatch From 6deedcc8d0fd3d299d859f5ec6a8d798188e7507 Mon Sep 17 00:00:00 2001 From: Anton Reinhard Date: Tue, 21 May 2024 10:15:49 +0200 Subject: [PATCH 18/22] Add (possibly temporary) fix for tests --- src/cross_sections.jl | 8 ++++---- src/phase_spaces.jl | 8 ++++++++ src/probabilities.jl | 8 ++++---- 3 files changed, 16 insertions(+), 8 deletions(-) diff --git a/src/cross_sections.jl b/src/cross_sections.jl index bf3279a..ab87c3e 100644 --- a/src/cross_sections.jl +++ b/src/cross_sections.jl @@ -143,8 +143,8 @@ function unsafe_differential_cross_section(phase_space_point::PhaseSpacePoint) phase_space_point.proc, phase_space_point.model, phase_space_point.ps_def, - momentum.(phase_space_point.in_particles), - momentum.(phase_space_point.out_particles), + momenta(phase_space_point, Incoming()), + momenta(phase_space_point, Outgoing()), ) end @@ -280,8 +280,8 @@ function differential_cross_section(phase_space_point::PhaseSpacePoint) phase_space_point.proc, phase_space_point.model, phase_space_point.ps_def, - momentum.(phase_space_point.in_particles), - momentum.(phase_space_point.out_particles), + momenta(phase_space_point, Incoming()), + momenta(phase_space_point, Outgoing()), ) end diff --git a/src/phase_spaces.jl b/src/phase_spaces.jl index 10590e2..724338b 100644 --- a/src/phase_spaces.jl +++ b/src/phase_spaces.jl @@ -321,6 +321,14 @@ struct PhaseSpacePoint{ end end +function momenta(psp::PhaseSpacePoint, ::Incoming) + return QEDbase._as_svec(momentum.(psp.in_particles)) +end + +function momenta(psp::PhaseSpacePoint, ::Outgoing) + return QEDbase._as_svec(momentum.(psp.out_particles)) +end + # recursion termination: success @inline _recursive_type_check(t::Tuple{}, p::Tuple{}, dir::ParticleDirection) = nothing diff --git a/src/probabilities.jl b/src/probabilities.jl index b5f24e0..19c2405 100644 --- a/src/probabilities.jl +++ b/src/probabilities.jl @@ -92,8 +92,8 @@ function unsafe_differential_probability(phase_space_point::PhaseSpacePoint) phase_space_point.proc, phase_space_point.model, phase_space_point.ps_def, - momentum.(phase_space_point.in_particles), - momentum.(phase_space_point.out_particles), + momenta(phase_space_point, Incoming()), + momenta(phase_space_point, Outgoing()), ) end @@ -230,8 +230,8 @@ function differential_probability(phase_space_point::PhaseSpacePoint) phase_space_point.proc, phase_space_point.model, phase_space_point.ps_def, - momentum.(phase_space_point.in_particles), - momentum.(phase_space_point.out_particles), + momenta(phase_space_point, Incoming()), + momenta(phase_space_point, Outgoing()), ) end From 3b479aba8b2e9b55f642654081a020bac40b9109 Mon Sep 17 00:00:00 2001 From: Anton Reinhard Date: Tue, 21 May 2024 10:33:37 +0200 Subject: [PATCH 19/22] Fix doctests --- src/phase_spaces.jl | 20 +++++++++----------- 1 file changed, 9 insertions(+), 11 deletions(-) diff --git a/src/phase_spaces.jl b/src/phase_spaces.jl index 724338b..0dd91d2 100644 --- a/src/phase_spaces.jl +++ b/src/phase_spaces.jl @@ -147,12 +147,10 @@ julia> using QEDbase; using QEDprocesses julia> ParticleStateful(Incoming(), Electron(), SFourMomentum(1, 0, 0, 0)) ParticleStateful: incoming electron - spin: all spins momentum: [1.0, 0.0, 0.0, 0.0] -julia> ParticleStateful(Outgoing(), Photon(), SFourMomentum(1, 0, 0, 0), PolX()) +julia> ParticleStateful(Outgoing(), Photon(), SFourMomentum(1, 0, 0, 0)) ParticleStateful: outgoing photon - polarization: x-polarized momentum: [1.0, 0.0, 0.0, 0.0] ``` """ @@ -225,25 +223,25 @@ julia> PhaseSpacePoint( Compton(), PerturbativeQED(), PhasespaceDefinition(SphericalCoordinateSystem(), ElectronRestFrame()), - [ + ( ParticleStateful(Incoming(), Electron(), SFourMomentum(1, 0, 0, 0)), ParticleStateful(Incoming(), Photon(), SFourMomentum(1, 0, 0, 0)) - ], - [ + ), + ( ParticleStateful(Outgoing(), Electron(), SFourMomentum(1, 0, 0, 0)), ParticleStateful(Outgoing(), Photon(), SFourMomentum(1, 0, 0, 0)) - ] + ) ) PhaseSpacePoint: process: one-photon Compton scattering model: perturbative QED phasespace definition: spherical coordinates in electron rest frame incoming particles: - -> incoming electron (all spins): [1.0, 0.0, 0.0, 0.0] - -> incoming photon (all polarizations): [1.0, 0.0, 0.0, 0.0] + -> incoming electron: [1.0, 0.0, 0.0, 0.0] + -> incoming photon: [1.0, 0.0, 0.0, 0.0] outgoing particles: - -> outgoing electron (all spins): [1.0, 0.0, 0.0, 0.0] - -> outgoing photon (all polarizations): [1.0, 0.0, 0.0, 0.0] + -> outgoing electron: [1.0, 0.0, 0.0, 0.0] + -> outgoing photon: [1.0, 0.0, 0.0, 0.0] ``` """ struct PhaseSpacePoint{ From a3ed26144057cb469dd781573e0cbe2effd59419 Mon Sep 17 00:00:00 2001 From: Anton Reinhard Date: Tue, 21 May 2024 10:57:21 +0200 Subject: [PATCH 20/22] Remove unused functions --- src/interfaces/process_interface.jl | 7 ------- 1 file changed, 7 deletions(-) diff --git a/src/interfaces/process_interface.jl b/src/interfaces/process_interface.jl index c889306..b790692 100644 --- a/src/interfaces/process_interface.jl +++ b/src/interfaces/process_interface.jl @@ -197,13 +197,6 @@ Return the number of outgoing particles of a given process. return length(outgoing_particles(proc_def)) end -@inline function incoming_particle_types(proc_def::AbstractProcessDefinition)::Type - return typeof(incoming_particles(proc_def)) -end -@inline function outgoing_particle_types(proc_def::AbstractProcessDefinition)::Type - return typeof(outgoing_particles(proc_def)) -end - """ in_phase_space_dimension( proc::AbstractProcessDefinition, From aaab7771f199b5db17bad644e53f95d901c2b112 Mon Sep 17 00:00:00 2001 From: Anton Reinhard Date: Tue, 21 May 2024 11:05:35 +0200 Subject: [PATCH 21/22] Make pretty print tests prettier --- test/phase_spaces.jl | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/test/phase_spaces.jl b/test/phase_spaces.jl index cd6d04f..40407fd 100644 --- a/test/phase_spaces.jl +++ b/test/phase_spaces.jl @@ -105,19 +105,19 @@ end @testset "Error handling" begin if (VERSION >= v"1.8") # julia versions before 1.8 did not have support for regex matching in @test_throws - @test_throws r"expected (.*)[Ii]ncoming(.*) (.*)[Pp]hoton(.*) but got (.*)[Oo]utgoing(.*) (.*)[Pp]hoton(.*)" PhaseSpacePoint( + @test_throws r"expected incoming photon but got outgoing photon" PhaseSpacePoint( process, model, phasespace_def, in_particles_invalid, out_particles_valid ) - @test_throws r"expected (.*)[Oo]utgoing(.*) (.*)[Pp]hoton(.*) but got (.*)[Ii]ncoming(.*) (.*)[Pp]hoton(.*)" PhaseSpacePoint( + @test_throws r"expected outgoing photon but got incoming photon" PhaseSpacePoint( process, model, phasespace_def, in_particles_valid, out_particles_invalid ) - @test_throws r"expected (.*)[Ii]ncoming(.*) (.*)[Ee]lectron(.*) but got (.*)[Ii]ncoming(.*) (.*)[Pp]hoton(.*)" PhaseSpacePoint( + @test_throws r"expected incoming electron but got incoming photon" PhaseSpacePoint( process, model, phasespace_def, (in_ph, in_el), out_particles_valid ) - @test_throws r"expected (.*)[Oo]utgoing(.*) (.*)[Ee]lectron(.*) but got (.*)[Oo]utgoing(.*) (.*)[Pp]hoton(.*)" PhaseSpacePoint( + @test_throws r"expected outgoing electron but got outgoing photon" PhaseSpacePoint( process, model, phasespace_def, in_particles_valid, (out_ph, out_el) ) end From 179fb0b6ac824800d1cc77b0c4773fff1235da05 Mon Sep 17 00:00:00 2001 From: Anton Reinhard Date: Tue, 21 May 2024 22:59:48 +0200 Subject: [PATCH 22/22] Review changes --- src/QEDprocesses.jl | 1 - src/phase_spaces.jl | 6 +++--- 2 files changed, 3 insertions(+), 4 deletions(-) diff --git a/src/QEDprocesses.jl b/src/QEDprocesses.jl index 871f56a..5aaf9cf 100644 --- a/src/QEDprocesses.jl +++ b/src/QEDprocesses.jl @@ -36,7 +36,6 @@ export AbstractFrameOfReference, CenterOfMomentumFrame, ElectronRestFrame export AbstractPhasespaceDefinition, PhasespaceDefinition export ParticleStateful, PhaseSpacePoint export spin, polarization, particle_direction, particle_species, momentum, momenta, getindex -export PhaseSpacePoint # specific compute models export PerturbativeQED diff --git a/src/phase_spaces.jl b/src/phase_spaces.jl index 0dd91d2..5010c4f 100644 --- a/src/phase_spaces.jl +++ b/src/phase_spaces.jl @@ -186,15 +186,15 @@ particle_species(part::ParticleStateful) = part.species momentum(part::ParticleStateful) = part.mom # recursion termination: base case -@inline _assemble_tuple_type(proc::Tuple{}, dir::ParticleDirection, ELTYPE::Type) = () +@inline _assemble_tuple_type(::Tuple{}, ::ParticleDirection, ::Type) = () # function assembling the correct type information for the tuple of ParticleStatefuls in a phasespace point constructed from momenta @inline function _assemble_tuple_type( - proc::Tuple{SPECIES_T,Vararg{AbstractParticleType}}, dir::DIR_T, ELTYPE::Type + particle_types::Tuple{SPECIES_T,Vararg{AbstractParticleType}}, dir::DIR_T, ELTYPE::Type ) where {SPECIES_T<:AbstractParticleType,DIR_T<:ParticleDirection} return ( ParticleStateful{DIR_T,SPECIES_T,ELTYPE}, - _assemble_tuple_type(proc[2:end], dir, ELTYPE)..., + _assemble_tuple_type(particle_types[2:end], dir, ELTYPE)..., ) end