Tools for working with EPREM simulation runs
Project description
EPREM(py)
Tools for working with EPREM simulation runs
Installation
$ pip install eprempy
Usage
Before getting started in earnest with eprempy
, one simple task you can perform to make sure everything is installed is to print the version number.
import eprempy
print(f"Using version {eprempy.__version__} of eprempy {chr(0x1f60a)}")
Using version 0.1.4 of eprempy 😊
Creating an observer
The examples below assume a working knowledge of EPREM — essentially, they assume that you have successfully run EPREM and have access to both the runtime configuration file and the output files.
Most eprempy
applications will begin by creating one or more observers. EPREM stores output in the form of "stream observers", with each observer file representing a set of connected simulation nodes. An observer file name has the form obsNNNNNN.nc
or fluxNNNNNN.nc
, where NNNNNN
is the zero-padded number of the EPREM node-stream, and the letters obs
or flux
indicate whether the corresponding data contains the full particle distribution or pre-computed flux.
To create an observer, the user must provide the observer ID, which is simply the EPREM stream number in the case of a stream observer. There are two additional attributes that are necessary for uniquely specifying an observer: the path to the directory containing the EPREM output and the path to the EPREM configuration file. When working within the directory containing EPREM output and the configuration file, the eprem.Observer
class will automatically detect the appropriate observer file and will search for a configuration file with a name matching a known pattern. However, you may always specify the data directory via source
and the name of the configuration file via config
. Either path may be relative to the current directory and may contain the ~
alias; the config
path may also be relative to the source
path.
from eprempy import eprem
stream = eprem.stream(4, source='data/with-dist/')
stream
Observer(4)
Each observer knows the location of its data file, as well as the runtime configuration file used to create the data.
print(stream.dataset.source.name)
print(stream.config.source.name)
obs000004.nc
eprem.cfg
We can request the value of simulation runtime parameters by aliased keyword. For example, let's check the assumed mean free path at 1 au. There are three possible ways we can refer to this parameter: 'lambda0', 'lamo', or 'lam0'.
for key in ('lambda0', 'lamo', 'lam0'):
print(f"{key} = {stream[key]!r}")
lambda0 = Variable([1.], unit='au')
lamo = Variable([1.], unit='au')
lam0 = Variable([1.], unit='au')
We can also request an array-like observable quantity by one of its aliases.
for alias in ('Vr', 'vr', 'Ur', 'ur'):
print(f"{alias} = {stream[alias]!r}")
Vr = Quantity(unit='m s^-1', dimensions={'time', 'shell'})
vr = Quantity(unit='m s^-1', dimensions={'time', 'shell'})
Ur = Quantity(unit='m s^-1', dimensions={'time', 'shell'})
ur = Quantity(unit='m s^-1', dimensions={'time', 'shell'})
As show above, each observable quantity knows its metric unit and array dimensions. Dimension names have the same meaning as in the raw EPREM output.
It is possible to change an observable quantity's unit by calling its withunit
method:
print(stream['vr'].unit)
print(stream['vr'].withunit('au / day').unit)
print(stream['flux'].unit)
print(stream['flux'].withunit('1 / (cm^2 s sr MeV)').unit)
m s^-1
au d^-1
J^-1 s^-1 sr^-1 m^-2
MeV^-1 s^-1 sr^-1 cm^-2
An observable quantity supports subscription by standard index types or by "measurable" indices. An index is measurable if it contains one or more numeric values followed by a unit-like string.
This is best illustrated by a concrete example. Let's first pick a particular time in hours and radius in au
time = stream['time'].withunit('hour')
print(f"time:\n{time[27]}\n")
radius = stream['radius'].withunit('au')
print(f"radius:\n{radius[27, 57]}")
time:
[67.2],
unit='h',
dimensions={'time'}
radius:
[[1.27798576]],
unit='au',
dimensions={'time', 'shell'}
Next, we'll retrieve the observable quantity representing plasma density.
rho = stream['rho']
print(f"density at fixed time:\n{rho[27, 56:59]}\n")
print(f"density at fixed radius:\n{rho[26:29, 57]}\n")
print(f"density at fixed time and radius:\n{rho[27, 57]}")
density at fixed time:
[[20884873.80632916 5146864.60956498 21476487.84191666]],
unit='m^-3',
dimensions={'time', 'shell'}
density at fixed radius:
[[ 5521952.49424706]
[ 5146864.60956498]
[18721643.43447628]],
unit='m^-3',
dimensions={'time', 'shell'}
density at fixed time and radius:
[[5146864.60956498]],
unit='m^-3',
dimensions={'time', 'shell'}
Finally, we'll create a measured index with value equivalent to our chosen time point and one equivalent to our chosen radius point. When we use the measured indices to subscript the density, we get a value very close to the value of density subscripted with standard indices.
t0 = 67.2, 'hour'
r0 = 1.27798576, 'au'
print(f"density at (t0, r0):\n{rho[t0, r0]}\n")
density at (t0, r0):
[[5146864.35972748]],
unit='m^-3',
dimensions={'time', 'radius'}
Of course, most of the steps above were simply for the sake of demonstration. You can acheive the desired affect in a single line
print(rho[(67.2, 'hour'), (1.27798576, 'au')])
[[5146864.35972748]],
unit='m^-3',
dimensions={'time', 'radius'}
Each observable quantity has a plot
method.
vr = stream['vr']
r0 = 0.5, 'au'
vr[:, r0].plot('k')
This method is useful for quickly viewing a single quantity. More complex plots require explicit use of a plotting module such as matplotlib.pyplot
. The following example creates a slightly more complex version of the above plot, and demonstrates use of the format
method for metric units.
import matplotlib.pyplot as plt
for r in (0.1, 0.5, 1.0):
plt.plot(time, vr[:, (r, 'au')], label=f"r = {r} au")
plt.legend()
plt.xlabel(f"Time [{time.unit}]")
plt.ylabel(rf"$V_r$ [{vr.unit.format(style='tex')}]")
plt.tight_layout()
You can print all observable quantities and parameters available to an observer by calling the observer's which
method.
stream.which('observables')
'phiOffset'
'times' == 't' == 'time'
'shells' == 'shell'
'pitch-angle cosine' == 'pitch angles' == 'pitch-angle cosines' == 'pitch-angles' == 'mu' == 'pitch angle' == 'pitch-angle'
'm' == 'mass'
'q' == 'charge'
'egrid' == 'energy' == 'energies' == 'E'
'speed' == 'v' == 'vgrid'
'R' == 'r' == 'radius'
'T' == 'theta'
'P' == 'phi'
'Br' == 'br'
'Btheta' == 'Bt' == 'btheta' == 'bt'
'bphi' == 'Bphi' == 'Bp' == 'bp'
'vr' == 'Vr' == 'Ur' == 'ur'
'utheta' == 'Vtheta' == 'ut' == 'Ut' == 'Vt' == 'Utheta'
'Uphi' == 'Up' == 'uphi' == 'up' == 'Vphi' == 'Vp'
'rho' == 'Rho'
'Dist' == 'dist' == 'f'
'J(E)' == 'flux' == 'Flux' == 'j' == 'J' == 'j(E)'
'x' == 'X'
'y' == 'Y'
'Z' == 'z'
'b mag' == 'b_mag' == 'b' == 'bmag' == 'B' == '|B|' == '|b|'
'|u|' == 'u_mag' == 'U' == '|U|' == 'u mag' == 'u' == 'umag'
'upara' == 'Upara' == 'u_para'
'uperp' == 'Uperp' == 'u_perp'
'angle' == 'flow angle' == 'flow_angle'
'divU' == 'divu' == 'div(U)' == 'div_u' == 'div U' == 'div(u)' == 'div u'
'rigidity' == 'R_g' == 'Rg'
'mean free path' == 'mean_free_path' == 'mfp'
'acceleration_rate' == 'ar' == 'acceleration rate'
'energy_density' == 'energy density'
'average energy' == 'average_energy'
'fluence'
'integral flux' == 'integral_flux' == 'intflux'
stream.which('parameters')
'minimum_energy' == 'Emin' == 'minimum energy'
'reference energy' == 'energy0'
'reference radius' == 'r0'
'FailModeDump'
'adiabaticChangeAlg'
'adiabaticFocusAlg'
'aziSunStart'
'boundaryFunctAmplitude' == 'J0'
'boundaryFunctBeta' == 'beta'
'E0' == 'boundaryFunctEcutoff'
'boundaryFunctGamma' == 'gamma'
'boundaryFunctXi' == 'xi'
'boundaryFunctionInitDomain'
'charge'
'checkSeedPopulation'
'dsh_hel_min'
'dsh_min'
'dumpFreq'
'dumpOnAbort'
'eMax'
'eMin'
'epCalcStartTime'
'epEquilibriumCalcDuration'
'epremDomain'
'epremDomainOutputTime'
'fieldAligned'
'flowMag'
'fluxLimiter'
'focusingLimit'
'gammaEhigh'
'gammaElow'
'idealShock'
'idealShockFalloff'
'idealShockInitTime'
'idealShockJump'
'idealShockPhi'
'idealShockScaleLength'
'idealShockSharpness'
'idealShockSpeed'
'idealShockTheta'
'idealShockWidth'
'idw_p'
'kperxkpar' == 'kper/kpar' == 'kper / kpar' == 'kper_kpar'
'lambda0' == 'lam0' == 'lamo'
'mass'
'mfpRadialPower' == 'mfp_radial_power'
'mhdBAu'
'mhdDensityAu'
'minInjectionEnergy'
'numColumnsPerFace'
'numEnergySteps'
'numEpSteps'
'numMuSteps'
'numNodesPerStream'
'numObservers'
'numRowsPerFace'
'numSpecies'
'obsPhi'
'obsR'
'obsTheta'
'omegaSun'
'outputFloat'
'outputRestart'
'parallelFlow'
'pointObserverOutput'
'pointObserverOutputTime'
'preEruptionDuration'
'rScale'
'rigidity_power' == 'rigidityPower'
'saveRestartFile'
'seedFunctionTest'
'shockDetectPercent'
'shockInjectionFactor'
'shockSolver'
'simStartTime'
'simStopTime'
'streamFluxOutput'
'streamFluxOutputTime'
'subTimeCouple'
'tDel'
'unifiedOutput'
'unifiedOutputTime'
'unstructuredDomain'
'unstructuredDomainOutputTime'
'useAdiabaticChange'
'useAdiabaticFocus'
'useBoundaryFunction'
'useDrift'
'useEPBoundary'
'useParallelDiffusion'
'useShellDiffusion'
'useStochastic'
'warningsFile'
Contributing
Interested in contributing? Check out the contributing guidelines. Please note that this project is released with a Code of Conduct. By contributing to this project, you agree to abide by its terms.
License
eprempy
was created by Matt Young. It is licensed under the terms of the BSD 3-Clause license.
Credits
eprempy
was created with cookiecutter
and the py-pkgs-cookiecutter
template.
Project details
Release history Release notifications | RSS feed
Download files
Download the file for your platform. If you're not sure which to choose, learn more about installing packages.