# Sympy codes to compute the inner expansion of Laplace and Helmholtz double-layer potentials on high aspect ratio ellipse

The expansions below are used in the manuscript

C. Carvalho, A. D. Kim, L. Lewis, and Z. Moitier, _Quadrature by Parity Asymptotic eXpansions (QPAX) for scattering by high aspect ratio particles_.

### Zoïs Moitier, Camille Carvalho (2021)
            
_Karlsruhe Institute of Technology, Germany_ 

_University of California Merced, USA_


```python
import sympy as sy
from IPython.display import display
```


```python
order = 3  # order of the expansion
```


```python
ε, δ, h = sy.symbols("varepsilon delta h", real=True, positive=True)
s, t, X = sy.symbols("s t X", real=True)
μ_term = sy.symbols(f"mu0:{order}", cls=sy.Function)
```

## The Kernel $K^L(s, t; \varepsilon)$ for Laplace's double-layer potential


```python
KL = -(ε) / (2 * sy.pi * (1 + ε ** 2 + (1 - ε ** 2) * sy.cos(s + t)))
display(KL)
```


$\displaystyle - \frac{\varepsilon}{2 \pi \left(\varepsilon^{2} + \left(1 - \varepsilon^{2}\right) \cos{\left(s + t \right)} + 1\right)}$


## Expansion of the kernel $K^L(s, \pi - s + \varepsilon X; \varepsilon)$


```python
KL_series = KL.subs(t, sy.pi - s + ε * X).series(x=ε, x0=0, n=order).expand()
display(KL_series)
```


$\displaystyle - \frac{1}{\pi X^{2} \varepsilon + 4 \pi \varepsilon} - \frac{X^{2} \varepsilon}{\pi X^{4} + 8 \pi X^{2} + 16 \pi} - \frac{X^{4} \varepsilon}{12 \pi X^{4} + 96 \pi X^{2} + 192 \pi} + O\left(\varepsilon^{3}\right)$


## Expansion of the function $\mu(\pi - s + \varepsilon X)$


```python
μ_series = sum((ε * X) ** i * μ_(sy.pi - s) for i, μ_ in enumerate(μ_term)) + sy.O(
    ε ** order
)
display(μ_series)
```


$\displaystyle \mu_{0}{\left(\pi - s \right)} + X \varepsilon \mu_{1}{\left(\pi - s \right)} + X^{2} \varepsilon^{2} \mu_{2}{\left(\pi - s \right)} + O\left(\varepsilon^{3}\right)$


## Expansion of the integrand $K^L(s, \pi-s+\varepsilon X; \varepsilon) \, \mu(\pi - s + \varepsilon X)\varepsilon$


```python
expr = (KL_series * μ_series * ε).expand()
display(expr)
```


$\displaystyle - \frac{\varepsilon \mu_{0}{\left(\pi - s \right)}}{\pi X^{2} \varepsilon + 4 \pi \varepsilon} - \frac{X \varepsilon^{2} \mu_{1}{\left(\pi - s \right)}}{\pi X^{2} \varepsilon + 4 \pi \varepsilon} - \frac{X^{2} \varepsilon^{2} \mu_{0}{\left(\pi - s \right)}}{\pi X^{4} + 8 \pi X^{2} + 16 \pi} - \frac{X^{2} \varepsilon^{3} \mu_{2}{\left(\pi - s \right)}}{\pi X^{2} \varepsilon + 4 \pi \varepsilon} - \frac{X^{4} \varepsilon^{2} \mu_{0}{\left(\pi - s \right)}}{12 \pi X^{4} + 96 \pi X^{2} + 192 \pi} + O\left(\varepsilon^{3}\right)$


## Inner expansion $\int_{-\frac{\delta}{2\varepsilon}}^{\frac{\delta}{2\varepsilon}} K^L(s, \pi-s+\varepsilon X; \varepsilon) \, \mu(\pi - s + \varepsilon X) \varepsilon\, \mathrm{d}X$


```python
int_expr = sy.O(ε ** 2)
for i in range(order):
    expr_term = expr.coeff(ε, i)
    for u_ in μ_term:
        expr_term_u = expr_term.coeff(u_(sy.pi - s), 1).factor()
        int_term = sy.integrate(expr_term_u, (X, -δ / (2 * ε), δ / (2 * ε))).expand()
        int_expr += (int_term * u_(sy.pi - s) * ε ** i).factor().expand()

display(int_expr)
```


$\displaystyle - \frac{\mu_{0}{\left(\pi - s \right)} \operatorname{atan}{\left(\frac{\delta}{4 \varepsilon} \right)}}{\pi} - \frac{\delta^{3} \varepsilon \mu_{0}{\left(\pi - s \right)}}{12 \pi \delta^{2} + 192 \pi \varepsilon^{2}} + O\left(\varepsilon^{2}\right)$


## The Kernel $K(s, t; \varepsilon) = K_1(s, t; \varepsilon)\log (4 \sin(\frac{s-t}{2})^2) + K_2(s, t; \varepsilon)$ for Helmholtz double-layer potential



```python
ε, δ, h, k = sy.symbols("varepsilon delta h k", real=True, positive=True)
s, t, X = sy.symbols("s t X", real=True)
u_term = sy.symbols(f"u0:{order}", cls=sy.Function)
```


```python
zε = (
    2
    * k
    * abs(sy.sin((s - t) / 2))
    * sy.sqrt(sy.cos((s + t) / 2) ** 2 + ε ** 2 * sy.sin((s + t) / 2) ** 2)
)
display(zε)
```


$\displaystyle 2 k \sqrt{\varepsilon^{2} \sin^{2}{\left(\frac{s}{2} + \frac{t}{2} \right)} + \cos^{2}{\left(\frac{s}{2} + \frac{t}{2} \right)}} \left|{\sin{\left(\frac{s}{2} - \frac{t}{2} \right)}}\right|$



```python
K1 = -(zε * KL * sy.besselj(1, zε)) / 2
display(K1)
```


$\displaystyle \frac{k \varepsilon \sqrt{\varepsilon^{2} \sin^{2}{\left(\frac{s}{2} + \frac{t}{2} \right)} + \cos^{2}{\left(\frac{s}{2} + \frac{t}{2} \right)}} \left|{\sin{\left(\frac{s}{2} - \frac{t}{2} \right)}}\right| J_{1}\left(2 k \sqrt{\varepsilon^{2} \sin^{2}{\left(\frac{s}{2} + \frac{t}{2} \right)} + \cos^{2}{\left(\frac{s}{2} + \frac{t}{2} \right)}} \left|{\sin{\left(\frac{s}{2} - \frac{t}{2} \right)}}\right|\right)}{2 \pi \left(\varepsilon^{2} + \left(1 - \varepsilon^{2}\right) \cos{\left(s + t \right)} + 1\right)}$



```python
log_term = sy.ln(4 * sy.sin((s - t) / 2) ** 2)
K2 = (
    -((sy.I * sy.pi * sy.hankel1(1, zε) + sy.besselj(1, zε) * log_term) * (KL * zε)) / 2
)
display(K2)
```


$\displaystyle \frac{k \varepsilon \sqrt{\varepsilon^{2} \sin^{2}{\left(\frac{s}{2} + \frac{t}{2} \right)} + \cos^{2}{\left(\frac{s}{2} + \frac{t}{2} \right)}} \left(\log{\left(4 \sin^{2}{\left(\frac{s}{2} - \frac{t}{2} \right)} \right)} J_{1}\left(2 k \sqrt{\varepsilon^{2} \sin^{2}{\left(\frac{s}{2} + \frac{t}{2} \right)} + \cos^{2}{\left(\frac{s}{2} + \frac{t}{2} \right)}} \left|{\sin{\left(\frac{s}{2} - \frac{t}{2} \right)}}\right|\right) + i \pi H^{(1)}_{1}\left(2 k \sqrt{\varepsilon^{2} \sin^{2}{\left(\frac{s}{2} + \frac{t}{2} \right)} + \cos^{2}{\left(\frac{s}{2} + \frac{t}{2} \right)}} \left|{\sin{\left(\frac{s}{2} - \frac{t}{2} \right)}}\right|\right)\right) \left|{\sin{\left(\frac{s}{2} - \frac{t}{2} \right)}}\right|}{2 \pi \left(\varepsilon^{2} + \left(1 - \varepsilon^{2}\right) \cos{\left(s + t \right)} + 1\right)}$


## Expansion of the kernel $K_1(s, \pi - s + \varepsilon X; \varepsilon)\log (4 \sin(\frac{s-t}{2})^2)u(s, \pi - s + \varepsilon X; \varepsilon)\varepsilon$



```python
def expansion(expr):
    return expr.subs(t, sy.pi - s + ε * X).series(x=ε, x0=0, n=order).expand()


dist = (
    sy.sqrt(2)
    * sy.Abs(sy.sin((s - t) / 2))
    * sy.sqrt(1 + ε ** 2 + (1 - ε ** 2) * sy.cos(s + t))
)
z_series = expansion(k * dist).subs(sy.cos(s) * sy.sign(sy.cos(s)), sy.Abs(sy.cos(s)))

J1 = sy.besselj(1, z_series.expand().removeO())
J1_series = J1.series(x=ε, x0=0, n=order)

u_series = sum((ε * X) ** i * u_(sy.pi - s) for i, u_ in enumerate(u_term)) + sy.O(
    ε ** order
)
log_term = sy.ln(1 + ε ** 2 + (1 - ε ** 2) * sy.cos(s + t))
log_series = expansion(log_term)

expr1 = (ε * z_series * u_series * J1_series * KL_series * log_series).expand()
display(expr1)
```


$\displaystyle - \frac{4 k^{2} \varepsilon^{3} \operatorname{u_{0}}{\left(\pi - s \right)} \log{\left(\varepsilon \right)} \cos^{2}{\left(s \right)}}{\pi X^{2} \varepsilon + 4 \pi \varepsilon} + \frac{4 k^{2} \varepsilon^{3} \operatorname{u_{0}}{\left(\pi - s \right)} \log{\left(2 \right)} \cos^{2}{\left(s \right)}}{2 \pi X^{2} \varepsilon + 8 \pi \varepsilon} - \frac{4 k^{2} \varepsilon^{3} \operatorname{u_{0}}{\left(\pi - s \right)} \log{\left(X^{2} + 4 \right)} \cos^{2}{\left(s \right)}}{2 \pi X^{2} \varepsilon + 8 \pi \varepsilon} - \frac{4 X k^{2} \varepsilon^{4} \operatorname{u_{1}}{\left(\pi - s \right)} \log{\left(\varepsilon \right)} \cos^{2}{\left(s \right)}}{\pi X^{2} \varepsilon + 4 \pi \varepsilon} - \frac{4 X k^{2} \varepsilon^{4} \operatorname{u_{0}}{\left(\pi - s \right)} \log{\left(\varepsilon \right)} \sin{\left(s \right)} \left|{\cos{\left(s \right)}}\right| \operatorname{sign}{\left(\cos{\left(s \right)} \right)}}{\pi X^{2} \varepsilon + 4 \pi \varepsilon} + \frac{4 X k^{2} \varepsilon^{4} \operatorname{u_{1}}{\left(\pi - s \right)} \log{\left(2 \right)} \cos^{2}{\left(s \right)}}{2 \pi X^{2} \varepsilon + 8 \pi \varepsilon} - \frac{4 X k^{2} \varepsilon^{4} \operatorname{u_{1}}{\left(\pi - s \right)} \log{\left(X^{2} + 4 \right)} \cos^{2}{\left(s \right)}}{2 \pi X^{2} \varepsilon + 8 \pi \varepsilon} + \frac{4 X k^{2} \varepsilon^{4} \operatorname{u_{0}}{\left(\pi - s \right)} \log{\left(2 \right)} \sin{\left(s \right)} \left|{\cos{\left(s \right)}}\right| \operatorname{sign}{\left(\cos{\left(s \right)} \right)}}{2 \pi X^{2} \varepsilon + 8 \pi \varepsilon} - \frac{4 X k^{2} \varepsilon^{4} \operatorname{u_{0}}{\left(\pi - s \right)} \log{\left(X^{2} + 4 \right)} \sin{\left(s \right)} \left|{\cos{\left(s \right)}}\right| \operatorname{sign}{\left(\cos{\left(s \right)} \right)}}{2 \pi X^{2} \varepsilon + 8 \pi \varepsilon} - \frac{X^{2} k^{2} \varepsilon^{3} \operatorname{u_{0}}{\left(\pi - s \right)} \log{\left(\varepsilon \right)} \cos^{2}{\left(s \right)}}{\pi X^{2} \varepsilon + 4 \pi \varepsilon} + \frac{X^{2} k^{2} \varepsilon^{3} \operatorname{u_{0}}{\left(\pi - s \right)} \log{\left(2 \right)} \cos^{2}{\left(s \right)}}{2 \pi X^{2} \varepsilon + 8 \pi \varepsilon} - \frac{X^{2} k^{2} \varepsilon^{3} \operatorname{u_{0}}{\left(\pi - s \right)} \log{\left(X^{2} + 4 \right)} \cos^{2}{\left(s \right)}}{2 \pi X^{2} \varepsilon + 8 \pi \varepsilon} - \frac{X^{3} k^{2} \varepsilon^{4} \operatorname{u_{1}}{\left(\pi - s \right)} \log{\left(\varepsilon \right)} \cos^{2}{\left(s \right)}}{\pi X^{2} \varepsilon + 4 \pi \varepsilon} - \frac{X^{3} k^{2} \varepsilon^{4} \operatorname{u_{0}}{\left(\pi - s \right)} \log{\left(\varepsilon \right)} \sin{\left(s \right)} \left|{\cos{\left(s \right)}}\right| \operatorname{sign}{\left(\cos{\left(s \right)} \right)}}{\pi X^{2} \varepsilon + 4 \pi \varepsilon} + \frac{X^{3} k^{2} \varepsilon^{4} \operatorname{u_{1}}{\left(\pi - s \right)} \log{\left(2 \right)} \cos^{2}{\left(s \right)}}{2 \pi X^{2} \varepsilon + 8 \pi \varepsilon} - \frac{X^{3} k^{2} \varepsilon^{4} \operatorname{u_{1}}{\left(\pi - s \right)} \log{\left(X^{2} + 4 \right)} \cos^{2}{\left(s \right)}}{2 \pi X^{2} \varepsilon + 8 \pi \varepsilon} + \frac{X^{3} k^{2} \varepsilon^{4} \operatorname{u_{0}}{\left(\pi - s \right)} \log{\left(2 \right)} \sin{\left(s \right)} \left|{\cos{\left(s \right)}}\right| \operatorname{sign}{\left(\cos{\left(s \right)} \right)}}{2 \pi X^{2} \varepsilon + 8 \pi \varepsilon} - \frac{X^{3} k^{2} \varepsilon^{4} \operatorname{u_{0}}{\left(\pi - s \right)} \log{\left(X^{2} + 4 \right)} \sin{\left(s \right)} \left|{\cos{\left(s \right)}}\right| \operatorname{sign}{\left(\cos{\left(s \right)} \right)}}{2 \pi X^{2} \varepsilon + 8 \pi \varepsilon} + O\left(\varepsilon^{4} \log{\left(\varepsilon \right)}\right)$


## Expansion of the kernel $K_2(s, \pi - s + \varepsilon X; \varepsilon)u(s, \pi - s + \varepsilon X; \varepsilon)\varepsilon$
In this case we expand the Hankel function of the first kind and first order using the NIST handbook of mathematical functions [https://dlmf.nist.gov/10.7]


```python
Y1approx = -2 / (sy.pi * z_series) + (2 / sy.pi) * sy.ln(z_series / 2) * J1_series
Y1_series = Y1approx.series(x=ε, x0=0, n=order)

expr2 = (
    ε
    * z_series
    * u_series
    * KL_series
    * (J1_series * log_series + sy.I * sy.pi * J1_series - sy.pi * Y1_series)
    / 2
).expand()
display(expr2)
```


$\displaystyle - \frac{\varepsilon \operatorname{u_{0}}{\left(\pi - s \right)}}{\pi X^{2} \varepsilon + 4 \pi \varepsilon} - \frac{X \varepsilon^{2} \operatorname{u_{1}}{\left(\pi - s \right)}}{\pi X^{2} \varepsilon + 4 \pi \varepsilon} + \frac{X \varepsilon^{2} \operatorname{u_{0}}{\left(\pi - s \right)} \sin{\left(s \right)} \left|{\cos{\left(s \right)}}\right| \operatorname{sign}{\left(\cos{\left(s \right)} \right)}}{2 \pi X^{2} \varepsilon \cos^{2}{\left(s \right)} + 8 \pi \varepsilon \cos^{2}{\left(s \right)}} - \frac{X \varepsilon^{2} \operatorname{u_{0}}{\left(\pi - s \right)} \sin{\left(s \right)} \operatorname{sign}{\left(\cos{\left(s \right)} \right)}}{2 \pi X^{2} \varepsilon \left|{\cos{\left(s \right)}}\right| + 8 \pi \varepsilon \left|{\cos{\left(s \right)}}\right|} + O\left(\varepsilon^{2}\right)$


## Inner expansion $\int_{-\frac{\delta}{2\varepsilon}}^{\frac{\delta}{2\varepsilon}} K_1(s, \pi-s+\varepsilon X; \varepsilon)\log (4 \sin(\frac{s-t}{2})^2) \, u(\pi - s + \varepsilon X) \varepsilon\, \mathrm{d}X$


```python
int_expr1 = sy.O(ε ** 2)
for i in range(order):
    expr1_term = expr1.coeff(ε, i)
    for u_ in u_term:
        expr1_term_u = expr1_term.coeff(u_(sy.pi - s), 1).factor()
        int_term1 = sy.integrate(expr1_term_u, (X, -δ / (2 * ε), δ / (2 * ε))).expand()
        int_expr1 += (int_term1 * u_(sy.pi - s) * ε ** i).factor().expand()

display(int_expr1)
```


$\displaystyle O\left(\varepsilon^{2}\right)$


## Inner expansion $\int_{-\frac{\delta}{2\varepsilon}}^{\frac{\delta}{2\varepsilon}} K_2(s, \pi-s+\varepsilon X; \varepsilon) u(\pi - s + \varepsilon X) \varepsilon\, \mathrm{d}X$



```python
int_expr2 = sy.O(ε ** 2)
for i in range(order):
    expr2_term = expr2.coeff(ε, i)
    for u_ in u_term:
        expr2_term_u = expr2_term.coeff(u_(sy.pi - s), 1).factor()
        int_term2 = sy.integrate(expr2_term_u, (X, -δ / (2 * ε), δ / (2 * ε))).expand()
        int_expr2 += (int_term2 * u_(sy.pi - s) * ε ** i).factor().expand()

display(int_expr2)
```


$\displaystyle - \frac{\operatorname{u_{0}}{\left(\pi - s \right)} \operatorname{atan}{\left(\frac{\delta}{4 \varepsilon} \right)}}{\pi} + O\left(\varepsilon^{2}\right)$


51dd244cbf712ff062faf548223e8c7ee8aa5bee
ipynb
Jupyter Notebook
Symbolic_expansions/SymPy/inner_expansions.ipynb
zmoitier/Scattering_BIE_QPAX
a4f9660f570c618f7a14585ab943b0f16c712632
MIT
Symbolic_expansions/SymPy/inner_expansions.ipynb
zmoitier/Scattering_BIE_QPAX
a4f9660f570c618f7a14585ab943b0f16c712632
MIT
Symbolic_expansions/SymPy/inner_expansions.ipynb
zmoitier/Scattering_BIE_QPAX
a4f9660f570c618f7a14585ab943b0f16c712632
MIT
Qwen/Qwen-72B
1. YES
2. YES