Feature/line integrator

biot_savart_integrators/generic/integrand.py

@@ -5,6 +5,7 @@
 """
 from scipy.constants import mu_0, pi
 from .utils import cross
+import xarray as xr
 
 
 _SI_FACTOR = mu_0 / (4*pi)
@@ -13,6 +14,16 @@
 def biot_savart_integrand(r, r_j, j, spatial_dim):
     R = r - r_j
     numerator = cross(j, R, spatial_dim)
-    denominator = (R**2).sum(dim=spatial_dim)**(3/2.)
+    denominator = (R**2).sum(dim=spatial_dim) ** (3 / 2.)
     integrand = _SI_FACTOR * numerator / denominator
-    return integrand
+    return integrand
+
+
+def biot_savart_potential_integrand(r: xr.DataArray, r_j: xr.DataArray,
+                                    j: xr.DataArray, spatial_dim: str):
+    R = r - r_j
+    numerator = j
+    denominator = (R ** 2).sum(dim=spatial_dim) ** (1 / 2.)
+    integrand = _SI_FACTOR * numerator / denominator
+    return integrand
+

biot_savart_integrators/generic/line_elements.py

@@ -5,11 +5,13 @@
 """
 import xarray as xr
 from .integrand import biot_savart_integrand as bsintegrand
+from .integrand import biot_savart_potential_integrand as bs_pot_integrands
+
 
 def biot_savart_integral(r: xr.DataArray, r_c: xr.DataArray,
                          dl: xr.DataArray, j: xr.DataArray,
-                         integration_dim: str, spatial_dim: str):
-    """Numerical integration of B-S law for set of line elements.
+                         integration_dim: str, spatial_dim: str) -> xr.DataArray:
+    """Numerical integration of B-S law for the set of line elements.
 
     :param r: Positions where the magnetic field is evaluated.
     :type r: xr.DataArray
@@ -34,3 +36,33 @@ def biot_savart_integral(r: xr.DataArray, r_c: xr.DataArray,
     # integral = integrand.sum(integration_dim)
     integral = integrand.integrate(integration_dim)
     return integral
+
+
+def biot_savart_potential_integral(r: xr.DataArray, r_c: xr.DataArray,
+                         dl: xr.DataArray, j: xr.DataArray,
+                         integration_dim: str, spatial_dim: str) -> xr.DataArray:
+    """Numerical integration of B-S law for vector potential for the set of line elements.
+
+    :param r: Positions where the magnetic field is evaluated.
+    :type r: xr.DataArray
+    :param r_c: Centers of line current elements.
+    :type r_c: xr.DataArray
+    :param dl: Lengt of the current element (scalar).
+    :type dl: xr.DataArray
+    :param j: Spatial current flowing through the element in Amps.
+    :type j: xr.DataArray
+    :param integration_dim: Dimension name over which index current
+                            elements.
+    :type integration_dim: str
+    :param spatial_dim: Name of the spetial dimension.
+    :type spatial_dim: str
+    """
+    j_int = j * dl
+    integrand = bs_pot_integrands(r, r_c, j_int, spatial_dim)
+
+    # TODO: Trapez rule used by integrate is more precise method then simple
+    # summation. However, this method may brake if the coordinates are assign
+    # to the integration_dim
+    # integral = integrand.sum(integration_dim)
+    integral = integrand.integrate(integration_dim)
+    return integral

tests/test_inf_wire.py

@@ -11,11 +11,12 @@
 @pytest.fixture
 def setup():
     I = 1.15766
-    a = 1.23
+    a = 1.0 # 1.23
     # high phi resolution to fulfill tolerance
 
     zs = xr.DataArray(np.linspace(-10000, 10000, int(1e6)), dims=['z'], name='z')
-    r0 = xr.DataArray([a, 0, 2.1], coords=[('s', list('xyz'))], name='r0')
+    #r0 = xr.DataArray([a, 0, 2.1], coords=[('s', list('xyz'))], name='r0')
+    r0 = xr.DataArray([a, 0, 0], coords=[('s', list('xyz'))], name='r0')
     r_c = xr.concat([xr.zeros_like(zs), xr.zeros_like(zs), zs], r0.s)
     return I, a, r0, r_c
 
@@ -35,3 +36,24 @@ def test_inf_wire_line_el(setup):
 
     np.testing.assert_allclose(B_line_el.sel(s=['x', 'z']), [0, 0])
     np.testing.assert_allclose(B_line_el.sel(s='y'), Btheta_analytic)
+
+
+def test_inf_wire_line_el_vector_potential(setup):
+    I, a, r0, r_c = setup
+    dl = r_c.differentiate("z")  # vector
+    j = dl * I  # vector with length element
+    dl = (dl ** 2).sum("s") ** 0.5  # actual lengths
+    j = j / dl  # normalization to line element length.
+
+    l = 10000
+    # Approximation for a << l
+    # https://phys.libretexts.org/Bookshelves/Electricity_and_Magnetism/Book%3A_Electricity_and_Magnetism_(Tatum)/09%3A_Magnetic_Potential/9.03%3A_Long%2C_Straight%2C_Current-carrying_Conductor
+
+    Az_analytical = (np.log(2*l) - np.log(a)) * mu_0 * I / (2 * np.pi)
+
+    A_calc = line_elements.biot_savart_potential_integral(r0, r_c, dl, j,
+                                                          integration_dim="z",
+                                                          spatial_dim="s")
+
+    np.testing.assert_allclose(A_calc.sel(s=['x', 'y']), [0, 0])
+    np.testing.assert_allclose(A_calc.sel(s='z'), Az_analytical)