In cosmic ray physics extensive progress has been made in recent years, both concerning theory and observation. The vast details in direct, indirect and secondary detections on the one hand provide the basis for a detailed modeling of the signatures via cosmic-ray transport and interaction, paving the way for the identification of Galactic cosmic-ray sources. On the other hand, the large number of constraints from different channels of cosmic-ray observables challenges these models frequently. In this review, we will summarize the state-of-the art of the detection of cosmic rays and their secondaries, followed by a discussion what we can learn from coupling our knowledge of the cosmic-ray observables to the theory of cosmic-ray transport in the Galactic magnetic field. Finally, information from neutral secondaries will be added to draw a multimessenger-picture of the non-thermal sky, in which the hypothesis of supernova remnants as the dominant sources survives best. While this has been known since the 1930s, evidence for this scenario is steadily growing, with the first detection of hadronic signatures at GeV energies detected for three SNRs with Fermi. The existence of SNRs as PeVatrons, however, is not validated yet. The discussion of this and other open questions concerning the level of anisotropy, composition and spectral shape of the cosmic-ray energy spectrum is reviewed. Future perspectives of how to find the smoking cosmic-ray source gun concludes this review.

For the first time, we characterize the rigidity regimes of the diffusion coefficient $\kappa$ for arbitrary rigidities and guide fields, which we derive as a function of physical and numerical parameters. We show that at turbulence levels above 5\% of the total magnetic field, the approximation of an energy dependence $\kappa\propto E^{1/3}$ as predicted for a Kolmogorov spectrum within Quasi-Linear Theory does not hold. Consequently, a proper description of cosmic-ray propagation can only be achieved by using a turbulence-level dependent diffusion coefficient and can contribute to solving the Galactic cosmic-ray gradient problem.

The popular JF12 analytic model by Jansson & Farrar provides a quantitative description of the Galaxy's large-scale magnetic field, which is widely used in various astrophysical applications. However, both the poloidal X-type component and the spiral disk component of JF12 exhibit regions in which the magnetic divergence constraint is violated. We first propose a cure for this problem, resulting in a truly solenoidal large-scale spiral field. Second, the otherwise straight field lines of the X-type component exhibit kinks in the Galactic plane that, in addition to implying the presence of a singular current sheet, may pose difficulties for e.g., numerical tracing of cosmic-ray particles. We propose and discuss two possible strategies to mitigate this problem. Although all corrections are kept as minimal as possible, the extended set of model parameters will have to be carefully readjusted in order to fully restore the agreement to observational data that the unmodified JF12 field is based on. Furthermore, the performance of our improved version of the field model is quantitatively assessed by test simulations using the CRPropa Galactic cosmic-ray propagation code.

The astrophysical interpretation of recent experimental observations of cosmic rays relies increasingly on Monte Carlo simulations of cosmic ray propagation and acceleration. Depending on the energy range of interest, several different propagation effects inside the Milky Way as well as in extragalactic space have to be taken into account when interpreting the data. With the CRPropa framework we aim to provide a toolbox for according simulations. In recent versions of CRPropa, the ballistic single particle propagation mode aiming primarily at extragalactic cosmic rays has been complemented by a solver for the differential transport equation to address propagation of galactic cosmic rays. Additionally, modules have been developed to address cosmic ray acceleration and many improvements have been added for simulations of electromagnetic secondaries. In this contribution we will give an overview of the CRPropa simulation framework with a focus on the latest improvements and highlight selected features by example applications.

Although several theories for the origin of cosmic rays in the region between the spectral `knee' and `ankle' exist, this problem is still unsolved. A variety of observations suggest that the transition from Galactic to extragalactic sources occurs in this energy range. In this work we examine whether a Galactic wind which eventually forms a termination shock far outside the Galactic plane can contribute as a possible source to the observed flux in the region of interest. Previous work by Bustard et al. (2017) estimated that particles can be accelerated up to energies above the `knee' up to $R_\mathrm{max} = 10^{16}$ eV for parameters drawn from a model of a Milky Way wind (Everett et al. 2017). A remaining question is whether the accelerated cosmic rays can propagate back into the Galaxy. To answer this crucial question, we simulate the propagation of the cosmic rays using the low energy extension of the CRPropa framework, based on the solution of the transport equation via stochastic differential equations. The setup includes all relevant processes, including three-dimensional anisotropic spatial diffusion, advection, and corresponding adiabatic cooling. We find that, assuming realistic parameters for the shock evolution, a possible Galactic termination shock can contribute significantly to the energy budget in the `knee' region and above. We estimate the resulting produced neutrino fluxes and find them to be below measurements from IceCube and limits by KM3NeT.

Radio galaxies are intensively discussed as the sources of cosmic rays observed above about $3\,{\times}\,10^{18}\,\text{eV}$, called ultra-high energy cosmic rays (UHECRs). We present a first, systematic approach that takes the individual characteristics of these sources into account, as well as the impact of the extragalactic magnetic-field structures up to a distance of 120 Mpc. We use a mixed simulation setup, based on 3D simulations of UHECRs ejected by observed, individual radio galaxies taken out to a distance of 120 Mpc, and on 1D simulations over a continuous source distribution contributing from beyond 120 Mpc. Additionally, we include the ultra-luminous radio galaxy Cygnus A at a distance of about $250\,$Mpc, as its contribution is so strong that it must be considered as an individual point source. The implementation of the UHECR ejection in our simulation setup is based on a detailed consideration of the physics of radio jets and standard first-order Fermi acceleration. We show that the average contribution of radio galaxies taken over a very large volume cannot explain the observed features of UHECRs measured at Earth. However, we obtain excellent agreement with the spectrum, composition, and arrival-direction distribution of UHECRs measured by the Pierre Auger Observatory, if we assume that most UHECRs observed arise from only two sources: The ultra-luminous radio galaxy Cygnus A, providing a mostly light composition of nuclear species dominating up to about $6\,{\times}\,10^{19}\,$eV, and the nearest radio galaxy Centaurus A, providing a heavy composition dominating above $6\,{\times}\,10^{19}\,$eV. Here we have to assume that extragalactic magnetic fields out to 250 Mpc, which we did not include in the simulation, are able to isotropize the UHECR events at about 8 EeV arriving from Cygnus A.

The propagation of charged cosmic rays through the Galactic environment influences all aspects of the observation at Earth. Energy spectrum, composition and arrival directions are changed due to deflections in magnetic fields and interactions with the interstellar medium. Today the transport is simulated with different simulation methods either based on the solution of a transport equation (multi-particle picture) or a solution of an equation of motion (single-particle picture). We developed a new module for the publicly available propagation software CRPropa 3.1, where we implemented an algorithm to solve the transport equation using stochastic differential equations. This technique allows us to use a diffusion tensor which is anisotropic with respect to an arbitrary magnetic background field. The source code of CRPropa is written in C++ with python steering via SWIG which makes it easy to use and computationally fast. In this paper, we present the new low-energy propagation code together with validation procedures that are developed to proof the accuracy of the new implementation. Furthermore, we show first examples of the cosmic ray density evolution, which depends strongly on the ratio of the parallel $\kappa_\parallel$ and perpendicular $\kappa_\perp$ diffusion coefficients. This dependency is systematically examined as well the influence of the particle rigidity on the diffusion process.

The luminosity ratio of electrons to protons as it is produced in stochastic acceleration processes in cosmic ray sources is an important quantity relevant for several aspects of the modeling of the sources themselves. It is usually assumed to be around 1:100 in the case of Galactic sources, while a value of 1:10 is typically assumed when describing extragalactic sources. It is supported by observations that the average ratios should be close to these values. At this point, however, there is no possibility to investigate how each individual source behaves. When looking at the physics aspects, a 1:100 ratio is well supported in theory when making the following assumptions: (1) the total number of electrons and protons that is accelerated are the same; (2) the spectral index of both populations after acceleration is $a_e=a_p\approx 2.2$. In this paper, we reinvestigate these assumptions. In particular, assumption (2) is not supported by observational data of the sources and PIC simulation yield different spectral indices as well. We present the detailed calculation of the electron-to-proton ratio, dropping the assumption of equal spectral indices. We distinguish between the ratio of luminosities and the ratio of the differential spectral behavior, which becomes necessary for cases where the spectral indices of the two particle populations are not the same. We discuss the possible range of values when allowing for different spectral indices concerning the spectral behavior of electrons and protons. Additionally, it is shown that the minimum energy of the accelerated population can have a large influence on the results.