Example of GT in the Magnetosphere

Consider the task of the calculation of the cut-off rigidity of antiprotons at altitude of 400 km, null meridian, and 40 degrees northern latitude in a vertical direction. These coordinates correspond to \(\vec{r} = (0.8139, 0, 0.6829)\) radius of earth, and the velocity along \(-\vec{r}\). Calculations are done in IGRF13 field. These initial parameters correspond to

Bfield = Gauss(model = "IGRF", model_type = core, version = 13)
InitialFlux = Flux(Spectrum = Monolines(energy = np.linspace(4.30, 6.00, 171) * Units.GeV),
                   Distribution = SphereSurf(Center = np.array([0.8139, 0, 0.6829]) * Units.RE, V0 = np.array([-0.8139, 0, -0.6829]), Radius=0)
                   Names = "apr",
                   Nevents = 171)

Overall the codes looks like:

import os
import numpy as np
from datetime import datetime

from gtsimulation.Global import Regions, Units
from gtsimulation.Algos import BunemanBorisSimulator

Region = Regions.Magnetosphere
Bfield = Gauss(model = "IGRF", model_type = core, version = 13)
Date = datetime(2006, 6, 15)

Medium = None

InitialFlux = Flux(Spectrum = Monolines(energy = np.linspace(4.30, 6.00, 171) * Units.GeV),
                   Distribution = SphereSurf(Center = np.array([0.8139, 0, 0.6829]) * Units.RE, V0 = np.array([-0.8139, 0, -0.6829]), Radius=0)
                   Names = "apr",
                   Nevents = 171)

UseDecay = False
NuclearInteraction = None

Nfiles = 1
Output = "IGRFSteps" + os.sep + "IGRF_Step5e-6"
Save = [4, {"Bfield":True}]
Verbose = True
BreakConditions = {"Rmin": 1 * Units.RE, "Rmax": 10 * Units.RE, "MaxPath": 1000 * Units.RE}
simulator = BunemanBorisSimulator(Date=Date,
                                  Region=Region,
                                  Bfield=Bfield,
                                  Particles=InitialFlux,
                                  Num=int(5e7),
                                  Step=5e-6,
                                  Save=Save,
                                  Nfiles=Nfiles,
                                  Output=Output,
                                  Verbose=Verbose,
                                  BreakCondition=BreakConditions,
                                  ForwardTrck=-1)
simulator()

The methodology of calculation of cut-off is the following. We backtrace the particle, and see whether it falls onto the earth, gets captured or leaves the magnetosphere. To turn on backtracking regime we pass in BunemanBorisSimulator(ForwardTrack=-1). These conditions are written in BreakConditions respectively.

BreakConditions = {"Rmin": 1 * Units.RE, "MaxPath": 1000 * Units.RE, "Rmax": 10 * Units.RE}

The particles that statisfy the first or the second condition, then these particles cannot be revered and the trajectories are forbidden. On the plot they are showed in black. Otherwise, the particles and reach the altitude from outside the magnetosphere, hence those are allowed trajectories (white lines on the plot).