iupitermag 0.2.5

Python package with Rust bindings to calculate the magnetic field in Jupiter's magnetosphere.

iupitermag

Introduction

iupitermag is a Python package to model Jupiter's magnetic field. This includes the internal field represented by spherical harmonics and the current sheet field like CON2020. In addition to using named field models like JRM33, JRM09, CON2020, you may also use your own model parameters to define a custom field model. It is very fast because it uses Rust internally.

Jupiter's internal magnetic field intensity at the 1-bar "surface" using the JRM33 model.

Many other public codes do this or something similar,

• JupiterMag - A Python package that uses the libjupitermag C++ library (both written by Matt James).
• jovian_jrm09_internal, con_2020 - Codes written for MATLAB and IDL by Marrisa Vogt.
• con2020 - Python code by Gabby Provan.

More details on Jupiter magnetic field models and these codes can be found in this paper - Wilson, R.J., Vogt, M.F., Provan, G. et al. Internal and External Jovian Magnetic Fields: Community Code to Serve the Magnetospheres of the Outer Planets Community. Space Sci Rev 219, 15 (2023). https://doi.org/10.1007/s11214-023-00961-3

Installation

Installing using wheels on PyPI

Because of breaking changes in numpy 2.0-2.4, the minimum supported Python version is 3.11. If you are on Python versions >=3.11, you can install directly using the wheels hosted on PyPI. This does not require you to install the Rust toolchain.

pip install iupitermag

If you are using a package manager like uv, you can also use the following commands. These will use the wheels uploaded on PyPI unless you are on a Python version or platform for which the wheels are unavailable.

uv pip install iupitermag

To add to uv.lock and pyproject.toml in your projects, use -

uv add iupitermag

Installing from source using uv

git clone https://github.com/ysar/iupitermag.git --depth=1
cd iupitermag/iupitermag-py
uv pip install .

Installing from source using maturin

git clone https://github.com/ysar/iupitermag.git --depth=1
cd iupitermag/iupitermag-py
maturin develop --release

Usage

Calculating the internal and current sheet fields at a single point.

All positions should be in the IAU_JUPITER coordinate system.

import numpy as np
import iupitermag as im

internal_field = im.InternalField("JRM33")
currentsheet_field = im.CurrentSheetField("CON2020")

r = 10.
theta = 0.
phi = 0.
b_int_rtp = internal_field.calc_field(r, theta, phi)
b_ext_rtp = currentsheet_field.calc_field(r, theta, phi)

# Or, you can the cartesian form.  This calls the spherical version internally.
x = 10.
y = 5.
z = 1.
b_int_xyz = internal_field.calc_field_xyz(x, y, z)
b_ext_xyz = currentsheet_field.calc_field_xyz(x, y, z)

Calculating the internal and current sheet fields for a collection of points.

If you have a collection of points stored as a single numpy array of shape (N, 3), you can use map_calc_field or parmap_calc_field (or their corresponding cartesian versions map_calc_field_xyz and parmap_calc_field_xyz).

# By default, iupitermag uses spherical coordinates (r, theta, phi)
points = np.zeros((10000, 3))
points[:, 0] = np.random.random_sample((10000,)) * 10 + 5

b_int = internal_field.map_calc_field(points)
b_ext = currentsheet_field.map_calc_field(points)

b_int = internal_field.parmap_calc_field(points)
b_ext = currentsheet_field.parmap_calc_field(points)

# Assume some cartesian coordinates in shape (N, 3)
points_xyz = points * 1.

b_int_xyz = internal_field.map_calc_field_xyz(points_xyz)
b_ext_xyz = currentsheet_field.map_calc_field_xyz(point_xyz)

b_int_xyz = internal_field.parmap_calc_field_xyz(points_xyz)
b_ext_xyz = currentsheet_field.parmap_calc_field_xyz(points_xyz)

The below figure shows the results of benchmarking different methods to calculate the JRM33 field for different number of points. Show are 1) a pure Python implementation of the internal field, 2) a Python loop that calls iupitermag's calc_field repeatedly on a number of points, 3) using map_calc_field directly on a points collections, and 4) using parmap_calc_field on a points collection, which is similar to map but uses Rayon for parallelization.

As you can see, parmap_ is usually the better option if you have more than 20 or so points, at least for this test that uses JRM33. Compared to the pure-Python implementation, parmap_ is about three orders of magnitude faster.

Tracing magnetic field lines

iupitermag can trace magnetic field lines to Jupiter using trace_field_to_planet, which takes as input a collection of starting points (which each result in a separate trace). The coordinates for these starting points should be cartesian.

import numpy as np
import iupitermag as im

internal_field = im.InternalField("JRM33")
currentsheet_field = im.CurrentSheetField("CON2020")

starting_positions_xyz = np.array([
    [-10., 0., 0.],
    [-15., 0., 0.],
    [-20., 0., 0.],
    [-25., 0., 0.],
    [10., 0., 0.],
    [15., 0., 0.],
    [20., 0., 0.],
    [25., 0., 0.],
])

trace = im.trace_field_to_planet(starting_positions_xyz, internal_field, currentsheet_field)

Using a custom internal and current sheet field

To use your own internal field or current sheet model, use the "Custom" argument when instantiating a field object. The following example describes how to create an internal field for an axially-aligned dipole and defining new current sheet model that is similar to CON2020 but ignores the current sheet tilt.

import numpy as np
import iupitermag as im

# Create a new axially aligned dipole field using Schmidt coefficients. 
# That is, only g[1, 0] is non-zero.
internal_field = im.InternalField(
    "Custom",
    g=np.array([[0.0, 0.0], [410993.4, 0.0]]),
    h=np.array([[0.0, 0.0], [0.0, 0.0]]),
)

# Get parameters for CON2020
con2020_field = im.CurrentSheetField("CON2020")
params = con2020_field.get_params()

# Change current sheet tilt to 0.0
params["theta_d"] = 0.0

# Define a modified current sheet field based on CON2020 parameters.
currentsheet_field = im.CurrentSheetField("Custom", params=params)

starting_positions = np.array(
    [
        [-10.0, 0.0, 0.0],
        [-15.0, 0.0, 0.0],
        [-20.0, 0.0, 0.0],
        [-25.0, 0.0, 0.0],
        [10.0, 0.0, 0.0],
        [15.0, 0.0, 0.0],
        [20.0, 0.0, 0.0],
        [25.0, 0.0, 0.0],
    ]
)

trace = im.trace_field_to_planet(
    starting_positions, internal_field, currentsheet_field
)