Computation and process setup interface (#14)

With this PR, the concept of *computation setups* is introduced, and
general as well as functionality related to scattering processes is
implemented.

## Description of the problem

One of the main tasks of this package is the computation of quantities
like differential and total cross-section for given scattering processes
and a given set of parameters. Usually, these quantities depend on a
fixed set of initial parameters, and their value needs to be calculated
for large amounts of input data. For example, the differential
cross-section of a process depends on the generic scattering process,
the compute model, and all sorts of initial parameters, which might be
physical (e.g. energy scales) or technical (e.g. integrator settings).
Once initialized, the differential cross-section is still a function of
the external momenta of a scattering process. Therefore, for a given set
of input data, i.e. momenta, one needs to be able to compute the
respective value of the differential cross-section using the given
setting.

## Suggested solution

The initial parameters for a given quantity will be collected in a
*setup*, which, once initialized, describes the *initialized quantity*,
i.e. the quantity with fixed initial parameters. The actual computation
can be performed by calling a member function `compute` on the setup
object and the respective input arguments. The implementation described
below defines an interface, that unifies this approach and extends the
computation workflow by adding input verification and post-processing
steps.

## Implementation details
The root type for all setups is `AbstractSetup`, for which the following
interface functions are defined:

```Julia 
_input_validation(stp::AbstractSetup, input::Any)
_compute(stp::AbstractSetup, input::Any)
_post_computation(stp::AbstractSetup, input::Any, result::Any)
```

None of those functions is exported, but they need to be added for a
concrete implementation of `AbstractSetup`. For all functions, except
`_compute`, a generic fallback is implemented, which uses a default
implementation for the respective function. Based on the interface
functions, the actual compute function is implemented as

```Julia

compute(stp::AbstractSetup, input::Any)

```
where internally, the following steps are performed:

1. input validation by calling `_input_validation` 
2. actual computation by calling `_compute`
3. post-processing by calling `_post_computation`

Additionally, based on `AbstractSetup`, a specialized version of setups
related to the combination of scattering processes and compute models is
implemented: `AbstractProcessSetup<:AbstractSetup`, where the following
functions are added to the setup interface:

```Julia
scattering_process(stp::AbstractProcessSetup)
compute_model(stp::AbstractProcessSetup)
```

Based on them, the following functions are delegated to the respective
type parameter:

```Julia
number_incoming_particles(stp::AbstractProcessSetup)
number_outgoing_particles(stp::AbstractProcessSetup)
```

## Final remarks

The possibilities for functionality based on the setup definition above
are not exhausted by the set of interface functions defined with this
PR. For example, the pre-and post-processing could be enhanced by more
fine granular structures. Furthermore, the initialization step of a
setup could have its input validation. This is currently delegated to
the users, who implement a setup by themself. Finally, the set of
functions defined on `AbstractProcessSetup` and delegated to the process
and model, could be extended if necessary. However, all the points above
should be considered only if they are necessary, which implies opening
dedicated issues if any.

---------

Co-authored-by: Uwe Hernandez Acosta <u.hernandez@hzdr.de>
Co-authored-by: Tom Jungnickel <140055258+tjungni@users.noreply.github.com>
Co-authored-by: Anton Reinhard <s1509337@msx.tu-dresden.de>

src/QEDprocesses.jl

@@ -13,6 +13,10 @@ export initial_phasespace_dimension, final_phasespace_dimension
 export number_incoming_particles, number_outgoing_particles
 export differential_cross_section, total_cross_section
 
+# Abstract setup interface
+export AbstractComputationSetup, InvalidInputError, compute
+export AbstractProcessSetup, scattering_process, physical_model
+
 # particle types
 export AbstractParticleType
 export FermionLike, Fermion, AntiFermion, MajoranaFermion
@@ -26,5 +30,6 @@ include("utils.jl")
 include("interfaces/particle_interface.jl")
 include("interfaces/model_interface.jl")
 include("interfaces/process_interface.jl")
+include("interfaces/setup_interface.jl")
 include("particle_types.jl")
 end

src/interfaces/model_interface.jl

@@ -3,10 +3,6 @@
 #
 # In this file, we define the interface of working with compute models in
 # general.
-# 
-# This file is part of `QEDprocesses.jl` which is by itself part of the `QED.jl`
-# ecosystem.
-#
 ###############
 # root type for models
 """

src/interfaces/particle_interface.jl

@@ -3,10 +3,6 @@
 #
 # In this file, we define the interface for working with particles in a general
 # sense. 
-# 
-# This file is part of `QEDprocesses.jl` which is by itself part of the `QED.jl`
-# ecosystem.
-#
 ###############
 
 

src/interfaces/process_interface.jl

@@ -3,10 +3,6 @@
 #
 # In this file, we define the interface for working with scattering processes in
 # general.
-# 
-# This file is part of `QEDprocesses.jl` which is by itself part of the `QED.jl`
-# ecosystem.
-#
 ###############
 
 """
@@ -47,7 +43,7 @@ function outgoing_particles end
 
 Return the number of incoming particles of a given process. 
 """
-@inline function number_incoming_pariticles(proc_def::AbstractProcessDefinition)
+@inline function number_incoming_particles(proc_def::AbstractProcessDefinition)
     return length(incoming_particles(proc_def))
 end
 
@@ -57,7 +53,7 @@ end
 
 Return the number of outgoing particles of a given process. 
 """
-@inline function number_outgoing_pariticles(proc_def::AbstractProcessDefinition)
+@inline function number_outgoing_particles(proc_def::AbstractProcessDefinition)
     return length(outgoing_particles(proc_def))
 end
 
@@ -70,8 +66,8 @@ end
         final_phasespace::AbstractVector{T},
     ) where {T<:QEDbase.AbstractFourMomentum}
 
-Interface function for the combination of scattering processes and models. Return the differential cross section of a 
-given process and model for a passed initial and final phase space. The elements of the `AbstractVector` representing the phase spaces 
+Interface function for the combination of scattering processes and physical models. Return the differential cross section of a 
+given process and physical model for a passed initial and final phase space. The elements of the `AbstractVector` representing the phase spaces 
 are the momenta of the respective particles. The implementation of this function for a concrete process and model must not 
 check if the length of the passed phase spaces match the respective number of particles. 
 
@@ -120,14 +116,14 @@ function differential_cross_section(
     init_phasespace::Union{AbstractVector{T},AbstractMatrix{T}},
     final_phasespace::Union{AbstractVector{T},AbstractMatrix{T}},
 ) where {T<:QEDbase.AbstractFourMomentum}
-    size(init_phasespace, 1) == number_incoming_pariticles(proc_def) || throw(
+    size(init_phasespace, 1) == number_incoming_particles(proc_def) || throw(
         DimensionMismatch(
-            "The number of momenta in the initial phasespace <{length(init_phasespace)}> does not match the number of incoming particles of the process <{number_incoming_pariticles(proc_def)}>.",
+            "The number of momenta in the initial phasespace <{length(init_phasespace)}> does not match the number of incoming particles of the process <{number_incoming_particles(proc_def)}>.",
         ),
     )
-    size(final_phasespace, 1) == number_outgoing_pariticles(proc_def) || throw(
+    size(final_phasespace, 1) == number_outgoing_particles(proc_def) || throw(
         DimensionMismatch(
-            "The number of momenta in the final phasespace <{length(final_phasespace)}> does not match the number of outgoing particles of the process <{number_outgoing_pariticles(proc_def)}>.",
+            "The number of momenta in the final phasespace <{length(final_phasespace)}> does not match the number of outgoing particles of the process <{number_outgoing_particles(proc_def)}>.",
         ),
     )
     return _differential_cross_section(
@@ -207,7 +203,7 @@ end
         init_phasespace::AbstractVector{T},
     ) where {T<:QEDbase.AbstractFourMomentum} end
 
-Interface function for the combination of scattering processes and models. Return the total cross section of a 
+Interface function for the combination of scattering processes and physical models. Return the total cross section of a 
 given process and model for a passed initial phase space. The elements of the `AbstractVector` representing the initial phase space
 are the momenta of the respective particles. The implementation of this function for a concrete process and model must not 
 check if the length of the passed initial phase spaces match number of incoming particles. 
@@ -255,7 +251,7 @@ end
         init_phasespace::Union{AbstractVector{T},AbstractMatrix{T}},
     ) where {T<:QEDbase.AbstractFourMomentum}
 
-Return the total cross section for a combination of a scattering process and a compute model evaluated on a given initial phase space. 
+Return the total cross section for a combination of a scattering process and a physical model evaluated on a given initial phase space. 
 
 This function will eventually call the respective interface function [`_total_cross_section`](@ref).
 
@@ -265,9 +261,9 @@ function total_cross_section(
     model_def::AbstractModelDefinition,
     init_phasespace::Union{AbstractVector{T},AbstractMatrix{T}},
 ) where {T<:QEDbase.AbstractFourMomentum}
-    size(init_phasespace, 1) == number_incoming_pariticles(proc_def) || throw(
+    size(init_phasespace, 1) == number_incoming_particles(proc_def) || throw(
         DimensionMismatch(
-            "The number of momenta in the initial phasespace <{length(init_phasespace)}> does not match the number of incoming particles of the process <{number_incoming_pariticles(proc_def)}>.",
+            "The number of momenta in the initial phasespace <{length(init_phasespace)}> does not match the number of incoming particles of the process <{number_incoming_particles(proc_def)}>.",
         ),
     )
     return _total_cross_section(proc_def, model_def, init_phasespace)

src/interfaces/setup_interface.jl

@@ -0,0 +1,194 @@
+###############
+# The process setup 
+#
+# In this file, we define the interface for general computation and process setups.
+###############
+
+"""
+Abstract base type for computation setups.  A *setup* means
+a collection of setup data needed to evaluate a dedicated quantity of given
+running data. Therefore, each setup is associated with a single quantity, which one may compute using the setup data and the running data. 
+Despite that, the decomposition into setup and running data is
+arbitrary, and this can be used for cases where a subset of the variables a
+quantity depends on is kept constant. 
+
+!!! note "Computation setup interface"
+    
+    The computation performed using a computation setup is separated into three steps:
+        
+        1. input validation
+        2. actual computation
+        3. post processing
+
+    where every step has its own interface function (see [`compute`](@ref) for details). 
+    
+    ## Input validation
+
+    Every subtype of `AbstractComputationSetup` should implement the interface function
+    
+    ```Julia
+    _assert_valid_input(stp::AbstractComputationSetup, input)
+    ```
+    
+    which should throw and an exception subtyped from [`AbstractInvalidInputException`](@ref) if the `input` is not valid for the computation of the associated quantity (see [`_assert_valid_input`](@ref) for more details). 
+    The default implementation does nothing, i.e. every input is valid by default. Provide a custom implementation if a different behavior is required.
+
+    ## Actual computation
+    
+    Every subtype of `AbstractComputationSetup` must at least implement the required interface function
+    
+    ```Julia
+    _compute(stp::AbstractComputationSetup, input) 
+    ```
+    
+    which computes the value of the associated quantity for a given `input` (see [`_compute`](@ref) for more details).
+
+
+    ## Post processing 
+
+    Every subtype of `AbstractComputationSetup` should implement the interface function
+    
+    ```Julia
+    _post_processing(stp::AbstractComputationSetup, input, result) 
+    ```
+    
+    which performs task after the actual computation, e.g. conversions or normalizations (see [`_post_processing`](@ref) for more details).
+
+"""
+abstract type AbstractComputationSetup end
+
+# convenience function to check if an object is a computation setup
+_is_computation_setup(::AbstractComputationSetup) = true
+
+"""
+Abstract base type for exceptions indicating invalid input. See [`InvalidInputError`](@ref) for a simple concrete implementation. 
+Concrete implementations should at least implement 
+
+```Julia
+
+Base.showerror(io::IO, err::CustomInvalidError) where {CustomInvalidError<:AbstractInvalidInputException}
+
+```
+"""
+abstract type AbstractInvalidInputException <: Exception end
+
+"""
+
+    InvalidInputError(msg::String)
+
+Exception which is thrown if a given input is invalid, e.g. passed to [`_assert_valid_input`](@ref).
+"""
+struct InvalidInputError <: AbstractInvalidInputException
+    msg::String
+end
+Base.showerror(io::IO, err::InvalidInputError) =
+    println(io, "InvalidInputError: $(err.msg).")
+
+"""
+
+    _assert_valid_input(stp::AbstractComputationSetup, input::Any)
+
+Interface function, which asserts that the given `input` is valid, and throws an [`InvalidInputError`](@ref) if not.
+
+!!! note "default implementation"
+
+    By default, every input is assumed to be valid. Therefore, this function does nothing. 
+    To customize this behavior, add your own implementation of
+
+    ```Julia
+    _assert_valid_input(stp::YourCustomSetup,input)
+    ```
+    which should throw an exception, which is a subtype of [`AbstractInvalidInputError`](@ref). One may also use the concrete implementation [`InvalidInputError`](@ref) if the input is invalid instead of writing a custom exception type.
+
+"""
+@inline function _assert_valid_input(stp::AbstractComputationSetup, input)
+    return nothing
+end
+
+"""
+    
+    function _post_processing(stp::AbstractComputationSetup, input::Any, result::Any)
+    
+Interface function, which is called in [`compute`](@ref) after [`_compute`](@ref) has been called. This function is dedicated to 
+finalize the result of a computation. 
+
+!!! note "default implementation"
+
+    Since in the case of no post processing the result of [`_compute`](@ref) is unchanged, this function returns `result` by default.
+
+"""
+@inline function _post_processing(stp::AbstractComputationSetup, input, result)
+    return result
+end
+
+"""
+    
+    _compute(stp::AbstractComputationSetup, input::Any)
+
+Interface function that returns the value of the associated quantity evaluated on `input`, which can be anything the associated quantity is defined to be feasible for.
+
+!!! note "unsafe implementation"
+
+    This function must be implemented for any subtype of [`AbstractComputationSetup`](@ref). It should not do any input validation or post processing (see [`_assert_valid_input`](@ref) and [`_post_processing`](@ref)), as those two are performed while calling 
+    the safe version of this function [`compute`](@ref).
+
+"""
+function _compute end
+
+"""
+
+    compute(stp::AbstractComputationSetup, input::Any)
+
+Return the value of the quantity associated with `stp` for a given `input`. 
+In addition to the actual call of the associated unsafe version [`_compute`](@ref),
+input validation ([`_assert_valid_input`]) and post processing 
+(using [`_post_processing`](@ref)) are wrapped around the calculation (see [`AbstractComputationSetup`](@ref) for details).
+"""
+function compute(stp::AbstractComputationSetup, input)
+    _assert_valid_input(stp, input)
+    raw_result = _compute(stp, input)
+    return _post_processing(stp, input, raw_result)
+end
+
+"""
+Abstract base type for setups related to combining scattering processes and physical models.  
+Every subtype of `AbstractProcessSetup` must implement at least the following 
+interface functions:
+
+```Julia
+scattering_process(::AbstractProcessSetup) 
+physical_model(::AbstractProcessSetup) 
+```
+
+Derived from these interface functions, the following delegations are provided:
+
+```Julia
+number_incoming_particles(::AbstractProcessSetup)
+number_outgoing_particles(::AbstractProcessSetup)
+```
+
+"""
+abstract type AbstractProcessSetup <: AbstractComputationSetup end
+
+"""
+    
+    scattering_process(stp::AbstractProcessSetup)
+
+Interface function that returns the scattering process associated with `stp`,
+i.e. an object which is a subtype of [`AbstractProcessDefinition`](@ref). 
+"""
+function scattering_process end
+
+"""
+
+    physical_model(stp::AbstractProcessSetup)
+
+Interface function that returns the physical model associated with `stp`, i.e.
+an object which is a subtype of [`AbstractModelDefinition`](@ref).
+"""
+function physical_model end
+
+@inline number_incoming_particles(stp::AbstractProcessSetup) =
+    number_incoming_particles(scattering_process(stp))
+@inline number_outgoing_particles(stp::AbstractProcessSetup) =
+    number_outgoing_particles(scattering_process(stp))

src/particle_types.jl

@@ -3,9 +3,6 @@
 #
 # In this file, we define the types of particles used in `QEDprocesses.jl` and
 # implement the abstact particle interface accordingly. 
-# 
-# This file is part of `QEDprocesses.jl` which is by itself part of the `QED.jl`
-# ecosystem.
 ###############
 
 """

src/utils.jl

@@ -2,9 +2,6 @@
 # utility functions
 #
 # This file contains small helper and utility functions used throughout the package.
-#
-# This file is part of `QEDprocesses.jl` which is by itself part of the `QED.jl`
-# ecosystem.
 ###############
 
 

test/interfaces/process_interface.jl

@@ -77,6 +77,11 @@ end
         @test outgoing_particles(TestProcess()) == OUTGOING_PARTICLES
     end
 
+    @testset "delegated functions" begin
+        @test number_incoming_particles(TestProcess()) == N_INCOMING
+        @test number_outgoing_particles(TestProcess()) == N_OUTGOING
+    end
+
 
     @testset "cross section" begin