Feature/line integrator

[kripnerl]

• add line_elements integrator
• add calculation of vector potential
• tests
  WIP: use integrate method vs sum:

[kripnerl]

This should be resolved before the merge. See the todo note.

[smartass101]

Ok, a simple solution would be to add an optional parameter integral_method='trapz' in a similar fashion as you do in pleque and add a warning about this into the docstring.

[smartass101]

I suppose that method='sum' would be something like intergand.sum(integration_dim) * dl?

[kripnerl]

Nope. integrate here only do "summation" of elements of the integrand, but with Trapez rule (see the example below). dl element is already included in j.

biot_savart_integrators/biot_savart_integrators/generic/curve.py

Line 12 in 32d2947

j = I*dl

If you assign some coordinates (for example l = cumsum(dl)) the integration will brake.

In [3]: a = xr.DataArray(np.linspace(1, 20, 20))                                                                                     

In [4]: a                                                                                                                            
Out[4]: 
<xarray.DataArray (dim_0: 20)>
array([ 1.,  2.,  3.,  4.,  5.,  6.,  7.,  8.,  9., 10., 11., 12., 13.,
       14., 15., 16., 17., 18., 19., 20.])
Dimensions without coordinates: dim_0

In [5]: a.sum()                                                                                                                      
Out[5]: 
<xarray.DataArray ()>
array(210.)

In [8]: a.integrate("dim_0")                                                                                                         
Out[8]: 
<xarray.DataArray ()>
array(199.5)

[kripnerl]

A possible solution would be something like this (untested; adjusted your function):

def biot_savart_integral(r, r_c, integration_dim, spatial_dim, I=1):
    # TODO make sure won't divide by coord - integration should cancel, but how accurately at edges?
    dl = r_c.differentiate(integration_dim)
    j = I*dl
    dl = (dl ** 2).sum(spatial_dim) ** 0.5 # normalize to be length element, not vector
    j = j / dl # remove length element

    integrand = bsintegrand(r, r_c, j, spatial_dim)
    integrand[integration_dim] = np.cumsum(dl)

    integral = integrand.integrate(integration_dim)
    return integral

or

def biot_savart_integral(r, r_c, integration_dim, spatial_dim, I=1):
    # TODO make sure won't divide by coord - integration should cancel, but how accurately at edges?
    dl = r_c.differentiate(integration_dim)
    j = I*dl
    integrand = bsintegrand(r, r_c, j, spatial_dim)
    
    integrand[integration_dim] = np.cumsum(dl)
    integral = integrand.integrate(integration_dim) / dl.sum()
    return integral

[smartass101]

Did you mean integrand.coords[integration_dim] = np.cumsum(dl) ?

I suppose the simplest solution would be just to use .sum() instead of .integrate(). Perhaps the possible inaccuracy of sum is not worth all these workarounds?

[smartass101]

Thank you for the proposed contribution. In general I think that .line_elements seems to be just a slight modification of .curve, so I think it would make more sense to merge it into one module with parametrized functions which could deliver both behaviors.

I'm sorry I didn't explicitly state that I prefer Numpydoc style docstrings, sorry. They are already present in some of the other modules in the ..torodial_sym package, so I would prefer to keep that convention.

[smartass101]

Ok, a simple solution would be to add an optional parameter integral_method='trapz' in a similar fashion as you do in pleque and add a warning about this into the docstring.

[kripnerl]

Ok I have added space in the notation. I am not an expert on the style. The type is in the declaration of the functions. But I can add typing to the docstring if you want.

[kripnerl]

There was quite major bug in _prepare_input_ function. It should ve fixed by 2fd95d1.

[kripnerl]

I've added there two more tests since I experienced some difficulties with signs. Everything is (and was) passing and can serve as and an example for calculation on some grid. I am not sure if this is the best way how to define grid in xarray though.