biot_savart

Follow Welling et al. (2020) calculate the magnetic field at Earth’s surface.

Integrate Biot-Savart over:

All currents within the Vlasov domain

Birkeland currents (FACs) in the “gap region” between the MHD inner boundary and the ionosphere, mapped (assuming J propto B) along the field lines connecting ionosphere radius R_IONO to coupling radius r_C

Horizontal Ionospheric currents (altitude 100 km)

The integration is a discrete summation over the 3 domains.

Integrating the domain #2 will likely require this mesh to be ‘refined’ in order to resolve the FAC structures see the functions refine_mesh() and graded_mesh()

This script has been tested with the Vlasiator runs EGL, FHA, FIA note that the usefulness for EGL is limited because there is no ionosphere (domain #3) for that run

To run this script, require access to the data reducer/variable ‘vg_J’ in the .vlsv file This may be supplied by a .vlsv file’s vlsvReader object, see keyword f_J_sidecar

This script is written for the UH environment. Adapt file paths as needed.

###

EXAMPLE CALL (1 thread): python biot_savart.py -nproc 1 -task 1 -run FHA

Example sidecar .vlsv files, containing ground magnetic field data, can be found at:: /wrk-vakka/group/spacephysics/vlasiator/3D/{run name}/sidecars/ig_B/

biot_savart.B_ionosphere(f, coord_list=None, ig_r=None, method='integrate')

Ionospheric (domain #3) current contribution to magnetic field

B_ionosphere() evaluates the magnetic field produced by ionospheric currents

Inputs:

f: VlsvReader object

keyword coord_list: locations where B-field is calculated

keyword ig_r: the ionospheric mesh

keyword method:

=’integrate’ (default): integrate over the whole ionosphere using Biot-Savart law

=’local’: ionosphere produces magnetic field locally, as by an infinite plane of current overhead

Outputs:

Magnetic field evaluated at coord_list

biot_savart.B_magnetosphere(f, f_J_sidecar=None, r_C=31855000.0, ig_r=None): Inner and outer magnetospheric contributions to Biot-Savart integral (domains #1, #2) wrapper for biot_savart()

biot_savart.b_dip(x, y, z, mag_mom_vector=array([0.e+00, 0.e+00, -8.e+22]))

Inputs: cartesian coordinates x,y,z [m]: keyword mag_mom_vector: Earth’s vector magnetic dipole moment

Outputs: dipole magnetic field

biot_savart.b_dip_direction(x, y, z, mag_mom_vector=array([0.e+00, 0.e+00, -8.e+22]))

Inputs: cartesian coordinates x,y,z [m]: keyword mag_mom_vector: Earth’s vector magnetic dipole moment

Outputs: dipole magnetic field unit vector

biot_savart.b_dip_magnitude(theta, r, mag_mom=8e+22)

Inputs: theta: colatitude [radians] r: radial distance [m] keyword mag_mom: magnetic dipole moment, default=8e22 [A / m^2], as in EGI, EGL, FIA, FHA runs

Outputs: Earth’s magnetic field magnitude [Tesla]

biot_savart.biot_savart(coord_list, f, f_J_sidecar=None, r_C=31855000.0, mesh='graded')

param coord_list: a list of 3-element arrays of coordinates [ [x1,y1,z1], [x2,y2,z2], … ], SI units: if considering just a single starting point, the code accepts a 3-element array-like object [x1,y1,z1]

f: vlsvReader object f_J_sidecar: vlsvReader object that contains pre-computed current ‘vg_J’

e.g., for EGL, files at: /wrk-vakka/group/spacephysics/vlasiator/3D/EGL/visualizations_2/ballooning/*.vlsv

runtime (FHA): overhead of about 200 sec (setup), plus 0.2 sec for each element of coord_list

returns a tuple (B_inner, B_outer) B_inner: the B-field [T] generated by FACs at 1 R_E < r< r_C B_outer: the B-field [T] generated by currents in simulation domain at r > r_C

biot_savart.cartesian_to_spherical(x, y, z)

r > 0 0 < theta < pi -pi < phi < pi all are assumed to be numpy arrays of equal dimensions

returns: r, theta, phi [tuple]

biot_savart.fac_map(f, vg_x, vg_y, vg_z, dx, f_J_sidecar=None, r_C=31855000.0, mag_mom_vector=array([0.e+00, 0.e+00, -8.e+22]))

Map the FACs along magnetic field lines (J propto B).

Inputs:
f: VlsvReader object f_J_sidecar: vlsvReader object that contains pre-computed current ‘vg_J’

if f_J_sidecar = None:
here is assumed the data reducer ‘ig_fac’ exists (such as for runs FHA and FIA) In this case, the currents in the FAC region (domain #2) will be mapped UP from the ionosphere ‘ig_’ grid

otherwise (f_J_sidecar = a *.vlsv string, for a file containing the current density J in the vlasov ‘vg_’ grid)
for run EGL, see files at /wrk-vakka/group/spacephysics/vlasiator/3D/EGL/visualizations_2/ballooning/*.vlsv In this case, the currents in the FAC region (domain #2) will be mapped DOWN from the vg_ grid

vg_x,vg_y,vg_z position [m], 1D numpy arrays.
note: these coordinates do not have to correspond with Vlasiator’s vg_ grid, This is relevant when fac_map() is called by biot_savart(), with keyword mesh=’refined’ or mesh=’graded’

dx is grid resolution [m], 1D numpy array, This is for the input cells which can be defined arbitrarily—not necessarily same as vg cells.

* coordinates, not just coordinates on the vg_ grid *

Returns:: the vector current density [A/m^2] at specified positions vg_x, vg_y, vg_z position output J=0 for cells that fall outside the FAC region (i.e if cell is not in region R_EARTH+dx/2 < r< r_C)

biot_savart.graded_mesh(x, y, z, dV, ns=array([8, 4, 2]), Rs=array([6371000., 12742000., 25484000.]))

Iteratively refine the mesh in ‘inner’ FAC region R_IONO < r < r_C, used to calculate FAC contribution to Biot-Savart integral

Calls refine mesh with different refinement factors n, at specified refinement boundaries ns: numpy array of refinement factors (ints) Rs: numpy array of refinement boundary radii [R_E] edge case: arrays are aligned so that for r>Rs[-1], use ns[-1] refinement

Want the refinement level to roughly scale as r^{-3/2}. To resolve the mapped features at all radii This works because magnetic dipole field goes as 1/r^3, roughly, and area A goes as r^2. (Need to set areas of the cell faces so B*A~constant, to resolve features mapped along field lines)

Example: graded_mesh(x, y, z, dV, ns = np.array([8, 4, 2]), Rs = np.array([R_EARTH, R_EARTH*2, R_EARTH*4])): –> refines 1D mesh size by factor of 8 for 1 RE<r<=2 RE, by factor of 4 for 2 RE<r<=4 RE, by factor of 2 for r>4 RE

biot_savart.integrate_biot_savart(coord_list, x, y, z, J, delta)

integration of the Biot-savart law magnetic field is evaluated at coordinates specified in coord_list (for example, the ionospheric coordinates)

Inputs:: x,y,z (1D array, size n): Cartesian coordinates delta (1D array, size n): the volume or area of the element being integrated J (2D array, size [3, n]): current density

all units SI
Outputs:: magnetic field evaluated at coord_list

Note: the units of J and delta depend on the type of integral (volume or surface), but the equation form is unchanged

Biot-Savart (volume): B = (mu_0 / 4 * pi) int { J x r’ / |r’|^3 } dV ([J] = A/m^2, delta == dV): (surface): B = (mu_0 / 4 * pi) int { J x r’ / |r’|^3 } dA ([J] = A/m, delta = dS)

biot_savart.mkdir_path(path): Make a directory from the stem of an input file name (path)

biot_savart.nearest_node_index(f, x, y, z, node_coords_iono=None): helper function for finding nearest ionospheric node

biot_savart.refine_mesh(x, y, z, dV, n)

refine the mesh in ‘inner’ FAC region R_IONO < r < r_C, used to calculate FAC contribution to Biot-Savart integral

x, y, z: initial coordinates (1D numpy float arrays), assumed to be at the center of cubic cell dV: cell volume 1D array or scalar value dx_0: side length of initial mesh (1D numpy array or single scalar value) n (int): the factor by which to refine the mesh

biot_savart.save_B_vlsv(input_tuple)

calculate magnetic field at the Earth’s surface and save in a .vslv file

Inputs: input tuple

input_tuple[0]: run (string) # ‘EGL’, ‘FHA’, or ‘FIA’ input_tuple[1]: fileIndex (int)

Outputs:

ig_r: Cartesian ionospheric grid locations (radius R_EARTH + 100km): note that B_iono, B_inner, B_outer are in fact evaluated at radius R_EARTH ig_r is a copy of the Vlasiator ionosphere mesh, to enable combination with other data reducers

B_iono: Ionospheric (Domain #3) contribution to ground magnetic field (radius R_EARTH) B_inner: Ionospheric (Domain #2) contribution to ground magnetic field (radius R_EARTH) B_outer: Ionospheric (Domain #1) contribution to ground magnetic field (radius R_EARTH)

Note: the input and output .vlsv file paths may need to be modified in this script for different users

biot_savart.spherical_to_cartesian(r, theta, phi)

r > 0 0 < theta < pi -pi < phi < pi all are assumed to be numpy arrays of equal dimensions

returns: x, y, z [tuple]

biot_savart.vec_len_2d(arr_2d): Vector length

biot_savart.vec_unit(arr_2d): assume arr_2d is a numpy array with shape [N, 3]. Return unit vectors with same shape