# Field Integration Techniques

Many analysis techniques for vector fields require solving an initial value problem for an arbitrary set of seed points and evaluating such solutions at a chosen resolution. Kamodo makes it easy to generate fieldline solutions by providing a function decorator that wraps scipy's powerful `solve_ivp` function. Each family of solutions is represented by a single function of a complex parameter. We illustrate the flexibility of this approach in the example below.


```python
# initialize
from plotly.offline import iplot, plot, init_notebook_mode
# init_notebook_mode(connected = True)

from kamodo import Kamodo, event, pointlike, kamodofy, solve
import numpy as np
import pandas as pd
```

## Dipole field model

We use the following dipole field model that can accept (m,) and (1,m), and (n,m) arrays.


```python
def Bdip(rvec):
    """Need math to work in a variety of arg shapes"""
    muvec = Bdip.muvec    
    r = np.linalg.norm(rvec, axis = 1)
    r[r==0] = np.nan

    try:
        rhat = rvec/r
    except:
        rhat = (rvec.T/r).T

    try:
        result = 3*np.dot(rhat, muvec.T)
    except:
        result = 3*np.dot(rhat.T, muvec.T).T


    result = (rhat.T*result).T

    try:
        result = result - muvec
    except:
        result = (result - muvec.T).T

    try:
        result = result/r**3
    except:
        result = (result.T/r**3).T

    return result

# set dipole moment
Bdip.muvec = np.array([0, 0, -1]) 

# pointlike enforces dimensionality
Bdip = pointlike(Bdip, '(n,m)->(n,m)', [np.float], squeeze = 0)
```


```python
kamodo = Kamodo()
```


```python
kamodo['Bvec'] = Bdip # register the dipole field
kamodo
```




\begin{equation}\vec{B}{\left(\vec{r} \right)} = \lambda{\left(\vec{r} \right)}\end{equation}



## Normalization

Instead of solving the initial value problem on the original field, we will be solving on the normalized field. This will mean that the integral path is the same as the arclength, allowing us to control the visual fidelity of the resulting field.

Create a normalization function to be applied to our field


```python
@kamodofy(equation = "$$\\hat{n}(\\vec{y}) = \\vec{y}/\\sqrt{\\vec{y} \\cdot \\vec{y}} $$")
@pointlike(signature = '(m,n)->(m,n)', squeeze = 0)
def normalized(yvec):   
    r = np.linalg.norm(yvec, axis = 1)
    r[r==0] = np.nan

    try:
        return yvec/r
    except:
        return (yvec.T/r).T


kamodo['nhat'] = normalized
```

Create a normalized field


```python
kamodo['bhat'] = "nhat(Bvec)"
kamodo
```




\begin{equation}\vec{B}{\left(\vec{r} \right)} = \lambda{\left(\vec{r} \right)}\end{equation}$$\hat{n}(\vec{y}) = \vec{y}/\sqrt{\vec{y} \cdot \vec{y}} $$\begin{equation}\hat{b}{\left(\vec{r} \right)} = \hat{n}{\left(\vec{B}{\left(\vec{r} \right)} \right)}\end{equation}



## Solving the initial value problem

Generate a set of seed points for integration


```python
x0 = np.linspace(-np.pi,np.pi,6)
y0 = np.linspace(-np.pi,np.pi,6)
z0 = 1

seeds = np.array(np.column_stack([c.ravel() for c in np.meshgrid(x0,y0,z0)]))
```

Create a stopping boundary for field line integrator


```python
@event
def boundary(s, rvec):
    r = np.linalg.norm(rvec)
    
    if np.isnan(r):
        result = 0
    else:
        result = r - 1
    return result

```

Solve the initial value problem for the normalized field


```python
kamodo['svec'] = solve(kamodo.bhat, # the field to be solved
                       seeds, # the initial positions
                       's', # the name of the integration parameter
                       (0,30), # the span to integrate over
                       npoints = 60, # the number of points to evaluate the solution
                       events = boundary, # stop at the boundary
                      )
kamodo
```




\begin{equation}\vec{B}{\left(\vec{r} \right)} = \lambda{\left(\vec{r} \right)}\end{equation}$$\hat{n}(\vec{y}) = \vec{y}/\sqrt{\vec{y} \cdot \vec{y}} $$\begin{equation}\hat{b}{\left(\vec{r} \right)} = \hat{n}{\left(\vec{B}{\left(\vec{r} \right)} \right)}\end{equation}\begin{equation}\vec{s}{\left(s \right)} = \lambda{\left(s \right)}\end{equation}



The solver returns a family of solutions, represented as a single function of a complex array, $\vec{s}(s)$ where
$s$ is a complex array. 

## Evaluating the Solutions
On evaluation, $\vec{s}(s)$ returns a pandas dataframe.


```python
kamodo.svec().head()
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th></th>
      <th>0</th>
      <th>1</th>
      <th>2</th>
    </tr>
    <tr>
      <th>seed</th>
      <th>integral</th>
      <th></th>
      <th></th>
      <th></th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th rowspan="5" valign="top">0.0</th>
      <th>-6.610169</th>
      <td>-0.347554</td>
      <td>-0.347554</td>
      <td>-0.924068</td>
    </tr>
    <tr>
      <th>-6.101695</th>
      <td>-0.615865</td>
      <td>-0.615865</td>
      <td>-1.261579</td>
    </tr>
    <tr>
      <th>-5.593220</th>
      <td>-0.922720</td>
      <td>-0.922720</td>
      <td>-1.525495</td>
    </tr>
    <tr>
      <th>-5.084746</th>
      <td>-1.256963</td>
      <td>-1.256963</td>
      <td>-1.713145</td>
    </tr>
    <tr>
      <th>-4.576271</th>
      <td>-1.608411</td>
      <td>-1.608411</td>
      <td>-1.822192</td>
    </tr>
  </tbody>
</table>
</div>



When using the default argument above, the solution evaluates at a resolution of npoints/span, stopping at the boundary.

## Complex parameterization
Kamodo represents the family of solutions to the initial value problem as a single function of a complex array.

The floor of the real part of the input parameter corresponds to the original seed array:


```python
kamodo.svec([0,1,2]).values
```




    array([[-3.14159265, -3.14159265,  1.        ],
           [-1.88495559, -3.14159265,  1.        ],
           [-0.62831853, -3.14159265,  1.        ]])



compare with original seeds:


```python
seeds[[0,1,2]]
```




    array([[-3.14159265, -3.14159265,  1.        ],
           [-1.88495559, -3.14159265,  1.        ],
           [-0.62831853, -3.14159265,  1.        ]])



The imaginary part denotes the integral along the corresponding solution. Here, we can choose evaluation points that were not in the original solution. Parameters outside the original span will be extrapolated.


```python
kamodo.svec([-6j, -5j, 0, 5j, 6j, 4 + 4j, 4 -5.777j])
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th></th>
      <th>0</th>
      <th>1</th>
      <th>2</th>
    </tr>
    <tr>
      <th>seed</th>
      <th>integral</th>
      <th></th>
      <th></th>
      <th></th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th rowspan="5" valign="top">0.0</th>
      <th>-6.000</th>
      <td>-0.674467</td>
      <td>-0.674467</td>
      <td>-1.320541</td>
    </tr>
    <tr>
      <th>-5.000</th>
      <td>-1.314574</td>
      <td>-1.314574</td>
      <td>-1.737228</td>
    </tr>
    <tr>
      <th>0.000</th>
      <td>-3.141593</td>
      <td>-3.141593</td>
      <td>1.000000</td>
    </tr>
    <tr>
      <th>5.000</th>
      <td>-0.120683</td>
      <td>-0.120683</td>
      <td>0.491788</td>
    </tr>
    <tr>
      <th>6.000</th>
      <td>0.115024</td>
      <td>0.115024</td>
      <td>-0.404619</td>
    </tr>
    <tr>
      <th rowspan="2" valign="top">4.0</th>
      <th>4.000</th>
      <td>0.094188</td>
      <td>-0.156979</td>
      <td>0.481816</td>
    </tr>
    <tr>
      <th>-5.777</th>
      <td>0.234822</td>
      <td>-0.391371</td>
      <td>-0.826931</td>
    </tr>
  </tbody>
</table>
</div>



## Plotting Fieldlines

We can quickly generate plots for all fieldlines at the default resolution by calling plot with the name of the fieldlines solution.


```python
import plotly.io as pio
```


```python
fig = kamodo.plot('svec')
pio.write_image(fig, 'images/fieldlines.svg')
```



To show the direction of the field at each point, we can evaluate $\hat{B}(\vec{s}(s))$


```python
fig = kamodo.plot('svec', 
                  Bhat = dict(rvec = kamodo.svec()))
pio.write_image(fig,'images/fieldlines_vectors.svg')
```



## Integration totals

To compute the total integral for each fieldline individually, we need a function to subtract the integration results at the endpoints.


```python
def integral(fieldline):
    endpoints = fieldline.reset_index().integral.iloc[[0,-1]]
    return endpoints.values[-1] - endpoints.values[0]
```


```python
totals = []
for seed, fieldline in kamodo.svec().groupby(level = 'seed'):
    totals.append(integral(fieldline))
    
totals[:5]
```




    [10.677966101694915,
     8.64406779661017,
     7.627118644067796,
     7.627118644067796,
     8.64406779661017]



Alternatively, we can use pandas' aggregation methods to apply our function on each fieldline.


```python
pd.DataFrame.groupby.a
```


```python
kamodo.svec().groupby(level='seed').aggregate(integral)
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>0</th>
      <th>1</th>
      <th>2</th>
    </tr>
    <tr>
      <th>seed</th>
      <th></th>
      <th></th>
      <th></th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>0.0</th>
      <td>10.677966</td>
      <td>10.677966</td>
      <td>10.677966</td>
    </tr>
    <tr>
      <th>1.0</th>
      <td>8.644068</td>
      <td>8.644068</td>
      <td>8.644068</td>
    </tr>
    <tr>
      <th>2.0</th>
      <td>7.627119</td>
      <td>7.627119</td>
      <td>7.627119</td>
    </tr>
    <tr>
      <th>3.0</th>
      <td>7.627119</td>
      <td>7.627119</td>
      <td>7.627119</td>
    </tr>
    <tr>
      <th>4.0</th>
      <td>8.644068</td>
      <td>8.644068</td>
      <td>8.644068</td>
    </tr>
    <tr>
      <th>5.0</th>
      <td>10.677966</td>
      <td>10.677966</td>
      <td>10.677966</td>
    </tr>
    <tr>
      <th>6.0</th>
      <td>8.644068</td>
      <td>8.644068</td>
      <td>8.644068</td>
    </tr>
    <tr>
      <th>7.0</th>
      <td>6.610169</td>
      <td>6.610169</td>
      <td>6.610169</td>
    </tr>
    <tr>
      <th>8.0</th>
      <td>5.084746</td>
      <td>5.084746</td>
      <td>5.084746</td>
    </tr>
    <tr>
      <th>9.0</th>
      <td>5.084746</td>
      <td>5.084746</td>
      <td>5.084746</td>
    </tr>
    <tr>
      <th>10.0</th>
      <td>6.610169</td>
      <td>6.610169</td>
      <td>6.610169</td>
    </tr>
    <tr>
      <th>11.0</th>
      <td>8.644068</td>
      <td>8.644068</td>
      <td>8.644068</td>
    </tr>
    <tr>
      <th>12.0</th>
      <td>7.627119</td>
      <td>7.627119</td>
      <td>7.627119</td>
    </tr>
    <tr>
      <th>13.0</th>
      <td>5.084746</td>
      <td>5.084746</td>
      <td>5.084746</td>
    </tr>
    <tr>
      <th>14.0</th>
      <td>5.593220</td>
      <td>5.593220</td>
      <td>5.593220</td>
    </tr>
    <tr>
      <th>15.0</th>
      <td>5.593220</td>
      <td>5.593220</td>
      <td>5.593220</td>
    </tr>
    <tr>
      <th>16.0</th>
      <td>5.084746</td>
      <td>5.084746</td>
      <td>5.084746</td>
    </tr>
    <tr>
      <th>17.0</th>
      <td>7.627119</td>
      <td>7.627119</td>
      <td>7.627119</td>
    </tr>
    <tr>
      <th>18.0</th>
      <td>7.627119</td>
      <td>7.627119</td>
      <td>7.627119</td>
    </tr>
    <tr>
      <th>19.0</th>
      <td>5.084746</td>
      <td>5.084746</td>
      <td>5.084746</td>
    </tr>
    <tr>
      <th>20.0</th>
      <td>5.593220</td>
      <td>5.593220</td>
      <td>5.593220</td>
    </tr>
    <tr>
      <th>21.0</th>
      <td>5.593220</td>
      <td>5.593220</td>
      <td>5.593220</td>
    </tr>
    <tr>
      <th>22.0</th>
      <td>5.084746</td>
      <td>5.084746</td>
      <td>5.084746</td>
    </tr>
    <tr>
      <th>23.0</th>
      <td>7.627119</td>
      <td>7.627119</td>
      <td>7.627119</td>
    </tr>
    <tr>
      <th>24.0</th>
      <td>8.644068</td>
      <td>8.644068</td>
      <td>8.644068</td>
    </tr>
    <tr>
      <th>25.0</th>
      <td>6.610169</td>
      <td>6.610169</td>
      <td>6.610169</td>
    </tr>
    <tr>
      <th>26.0</th>
      <td>5.084746</td>
      <td>5.084746</td>
      <td>5.084746</td>
    </tr>
    <tr>
      <th>27.0</th>
      <td>5.084746</td>
      <td>5.084746</td>
      <td>5.084746</td>
    </tr>
    <tr>
      <th>28.0</th>
      <td>6.610169</td>
      <td>6.610169</td>
      <td>6.610169</td>
    </tr>
    <tr>
      <th>29.0</th>
      <td>8.644068</td>
      <td>8.644068</td>
      <td>8.644068</td>
    </tr>
    <tr>
      <th>30.0</th>
      <td>10.677966</td>
      <td>10.677966</td>
      <td>10.677966</td>
    </tr>
    <tr>
      <th>31.0</th>
      <td>8.644068</td>
      <td>8.644068</td>
      <td>8.644068</td>
    </tr>
    <tr>
      <th>32.0</th>
      <td>7.627119</td>
      <td>7.627119</td>
      <td>7.627119</td>
    </tr>
    <tr>
      <th>33.0</th>
      <td>7.627119</td>
      <td>7.627119</td>
      <td>7.627119</td>
    </tr>
    <tr>
      <th>34.0</th>
      <td>8.644068</td>
      <td>8.644068</td>
      <td>8.644068</td>
    </tr>
    <tr>
      <th>35.0</th>
      <td>10.677966</td>
      <td>10.677966</td>
      <td>10.677966</td>
    </tr>
  </tbody>
</table>
</div>



9fd7a81dfe7aedc132494de966296e16b347184e
ipynb
Jupyter Notebook
docs/notebooks/FieldIntegration.ipynb
EnsembleGovServices/Kamodo-ccmc-readers
75841f7ad832997159046d4b2523e0a244316e9d
NASA-1.3
2019-08-30T16:18:53.000Z
2022-01-13T16:33:51.000Z
docs/notebooks/FieldIntegration.ipynb
EnsembleGovServices/Kamodo-ccmc-readers
75841f7ad832997159046d4b2523e0a244316e9d
NASA-1.3
2020-02-13T20:46:13.000Z
2022-01-18T17:00:15.000Z
docs/notebooks/FieldIntegration.ipynb
EnsembleGovServices/Kamodo-ccmc-readers
75841f7ad832997159046d4b2523e0a244316e9d
NASA-1.3
2019-08-16T21:22:40.000Z
2022-03-07T15:39:18.000Z
Qwen/Qwen-72B
1. YES
2. YES