#!/usr/bin/env python
# coding: utf-8

# ###### Content provided under a Creative Commons Attribution license CC-BY 4.0; code under BSD 3-Clause license. (c)2015 L.A. Barba, Pi-Yueh Chuang.

# # Exercise: Derivation of the vortex-source panel method

# The potential at location $(x, y)$ induced by an uniform flow, a source sheet, and a vortex sheet can be represented as

# $$
# \begin{equation}
# \begin{split}
# \phi(x, y) 
# &= \phi_{uniform\ flow}(x, y) \\ 
# &+ \phi_{source\ sheet}(x, y) + \phi_{vortex\ sheet}(x, y)
# \end{split}
# \end{equation}
# $$

# That is

# $$
# \begin{equation}
# \begin{split}
# \phi(x, y) &= xU_{\infty}\cos(\alpha) + yU_{\infty}\sin(\alpha) \\
# &+
# \frac{1}{2\pi} \int_{sheet} \sigma(s)\ln\left[(x-\xi(s))^2+(y-\eta(s))^2\right]^{\frac{1}{2}}ds \\
# &-
# \frac{1}{2\pi} \int_{sheet} \gamma(s)\tan^{-1} \frac{y-\eta(s)}{x-\xi(s)}ds
# \end{split}
# \end{equation}
# $$

# where $s$ is local coordinate on the sheet, and $\xi(s)$ and $\eta(s)$ are coordinate of the infinite source and vortex on the sheet. In the above equation, we assume the source sheet and the vortex sheet overlap.

# ------------------------------------------------------

# ### Q1:

# If we discretize the sheet into $N$ panels, re-write the above equation using discretized integral. Assume $l_j$ represents the length of the panel $j$. And so that
# 
# $$
# \begin{equation}
# \left\{
# \begin{array}{l}
# \xi_j(s)=x_j-s\sin\beta_j \\
# \eta_j(s)=y_j+s\cos\beta_j
# \end{array}
# ,\ \ \ 
# 0\le s \le l_j
# \right.
# \end{equation}
# $$
# 
# The following figure shows the panel $j$:
# 
# <center> <img src="resources/Lesson11_Exercise_Fig.1.png" width=360> </center>
# 
# HINT: for example, consider the integral $\int_0^L f(x) dx$, if we discretize the domain $0\sim L$ into 3 panels, the integral can be writen as:
# 
# $$
# \int_0^L f(x) dx = \int_0^{L/3} f(x)dx+\int_{L/3}^{2L/3} f(x)dx+\int_{2L/3}^{L} f(x)dx \\
# = \sum_{j=1}^3 \int_{l_j}f(x)dx
# $$

# ----------------------------

# Now let's assume
# 
# 1. $\sigma_j(s) = constant = \sigma_j$
# 2. $\gamma_1(s) = \gamma_2(s) = ... = \gamma_N(s) = \gamma$

# ------------------------------------------------

# ### Q2:

# Apply the above assumption into the equation of $\phi(x, y)$ you derived in Q1.

# ---------------------------

# The normal velocity $U_n$ can be derived from the chain rule:

# $$
# \begin{equation}
# \begin{split}
# U_n &= \frac{\partial \phi}{\partial \vec{n}} \\
# &=
# \frac{\partial \phi}{\partial x}\frac{\partial x}{\partial \vec{n}}
# +
# \frac{\partial \phi}{\partial y}\frac{\partial y}{\partial \vec{n}} \\
# &=
# \frac{\partial \phi}{\partial x}\nabla x\cdot \vec{n}
# +
# \frac{\partial \phi}{\partial y}\nabla y\cdot \vec{n} \\
# &=
# \frac{\partial \phi}{\partial x}n_x
# +
# \frac{\partial \phi}{\partial y}n_y
# \end{split}
# \end{equation}
# $$

# The tangential velocity can also be obtained using the same technique. So we can have the normal and tangential velocity at the point $(x, y)$ using:

# $$
# \begin{equation}
# \left\{
# \begin{array}{l}
# U_n(x, y)=\frac{\partial \phi}{\partial x}(x, y) n_x(x, y)+\frac{\partial \phi}{\partial y}(x, y) n_y(x, y) \\
# U_t(x, y)=\frac{\partial \phi}{\partial x}(x, y) t_x(x, y)+\frac{\partial \phi}{\partial y}(x, y) t_y(x, y)
# \end{array}
# \right.
# \end{equation}
# $$

# -------------------------------------

# ### Q3:

# Using the above equation, derive the $U_n(x,y)$ and $U_t(x,y)$ from the equation you obtained in Q2.

# -----------------------------------------

# ### Q4:

# Consider the normal velocity at the center of $i$-th panel, i.e., $(x_{c,i}, y_{c,i})$, after replacing $(x_{c,i}, y_{c,i})$ with $(x, y)$ in the equation you derived in the Q3, we can re-write the equation in matrix form: 

# $$
# \begin{equation}
# \begin{split}
# U_n(x_{c,i}, y_{c,i}) &= U_{n,i} \\
# &= b^n_i + \left[\begin{matrix} A^n_{i1} && A^n_{i2} && ... && A^n_{iN}\end{matrix}\right]\left[\begin{matrix} \sigma_1 \\ \sigma_2 \\ \vdots \\ \sigma_N \end{matrix}\right] + \left(\sum_{j=1}^N B^n_{ij}\right)\gamma \\
# &= b^n_i + \left[\begin{matrix} A^n_{i1} && A^n_{i2} && ... && A^n_{iN} && \left(\sum_{j=1}^N B^n_{ij}\right) \end{matrix}\right]\left[\begin{matrix} \sigma_1 \\ \sigma_2 \\ \vdots \\ \sigma_N \\ \gamma \end{matrix}\right]
# \end{split}
# \end{equation}
# $$

# $$
# \begin{equation}
# \begin{split}
# U_t(x_{c,i}, y_{c,i}) &= U_{t,i} \\
# &= b^t_i + \left[\begin{matrix} A^t_{i1} && A^t_{i2} && ... && A^t_{iN}\end{matrix}\right]\left[\begin{matrix} \sigma_1 \\ \sigma_2 \\ \vdots \\ \sigma_N \end{matrix}\right] + \left(\sum_{j=1}^N B^t_{ij}\right)\gamma \\
# &= b^t_i + \left[\begin{matrix} A^t_{i1} && A^t_{i2} && ... && A^t_{iN} && \left(\sum_{j=1}^N B^t_{ij}\right) \end{matrix}\right]\left[\begin{matrix} \sigma_1 \\ \sigma_2 \\ \vdots \\ \sigma_N \\ \gamma \end{matrix}\right]
# \end{split}
# \end{equation}
# $$

# What are the $b^n_i$, $A^n_{ij}$, $B^n_{ij}$, $b^t_i$, $A^t_{ij}$, and $B^t_{ij}$? 

# -----------------------

# Given the fact that (from the Fig. 1) 
# 
# $$
# \begin{equation}
# \left\{\begin{matrix} \vec{n}_i=n_{x,i}\vec{i}+n_{y,i}\vec{j} = \cos(\beta_i)\vec{i}+\sin(\beta_i)\vec{j} \\ \vec{t}_i=t_{x,i}\vec{i}+t_{y,i}\vec{j} = -\sin(\beta_i)\vec{i}+\cos(\beta_i)\vec{j} \end{matrix}\right.
# \end{equation}
# $$
# 
# we have
# 
# $$
# \begin{equation}
# \left\{
# \begin{matrix}
# n_{x,i}=t_{y,i} \\
# n_{y,i}=-t_{x,i}
# \end{matrix}
# \right.
# ,\ or\ 
# \left\{
# \begin{matrix}
# t_{x,i}=-n_{y,i} \\
# t_{y,i}=n_{x,i}
# \end{matrix}
# \right.
# \end{equation}
# $$

# -----------------------

# ### Q5:

# Applying the above relationship between $\vec{n}_i$ and $\vec{t}_i$ to your answer of the Q4, you should find that relationships exist between $B^n_{ij}$ and $A^t_{ij}$ and between $B^t_{ij}$ and $A^n_{ij}$. This means, in your codes, you don't have to actually calculate the $B^n_{ij}$ and $B^t_{ij}$. What are the relationship?

# -------------------------

# Now, note that when $i=j$, there is a singular point in the integration domain when calculating $A^n_{ii}$ and $A^t_{ii}$. This singular point occurs when $s=l_i/2$, i.e., $\xi_i(l_i/2)=x_{c,i}$ and $\eta_i(l_i/2)=y_{c,i}$. This means we need to calculate $A^n_{ii}$ and $A^t_{ii}$ analytically.

# --------------------------

# ### Q6:

# What is the exact values of $A^n_{ii}$ and $A^t_{ii}$?

# ------------------------------

# In our problem, there are $N+1$ unknowns, that is, $\sigma_1, \sigma_2, ..., \sigma_N, \gamma$. We'll need $N+1$ linear equations to solve the unknowns. The first $N$ linear equations can be obtained from the non-penetration condition on the center of each panel. That is
# 
# $$
# \begin{equation}
# \begin{split}
# U_{n,i} &= 0 \\
# &= b^n_i + \left[\begin{matrix} A^n_{i1} && A^n_{i2} && ... && A^n_{iN} && \left(\sum_{j=1}^N B^n_{ij}\right) \end{matrix}\right]\left[\begin{matrix} \sigma_1 \\ \sigma_2 \\ \vdots \\ \sigma_N \\ \gamma \end{matrix}\right] \\
# &,\ \ for\ i=1\sim N
# \end{split}
# \end{equation}
# $$
# 
# or
# 
# $$
# \begin{equation}
# \begin{split}
# &\left[\begin{matrix} A^n_{i1} && A^n_{i2} && ... && A^n_{iN} && \left(\sum_{j=1}^N B^n_{ij}\right) \end{matrix}\right]\left[\begin{matrix} \sigma_1 \\ \sigma_2 \\ \vdots \\ \sigma_N \\ \gamma \end{matrix}\right] =-b^n_i \\
# &,\ \ for\ i=1\sim N
# \end{split}
# \end{equation}
# $$

# For the last equation, we use Kutta-condition to obtain that.
# 
# $$
# \begin{equation}
# U_{t,1} = - U_{t,N}
# \end{equation}
# $$

# ----------------------

# ### Q7:

# Apply the matrix form of the $U_{t,i}$ and $U_{t,N}$ to the Kutta-condition and obtain the last linear equation. Re-arrange the equation so that unknowns are always on the LHS while the knowns on RHS.

# ---------------------

# ### Q8:

# Now you have $N+1$ linear equations and can solve the $N+1$ unknowns. Try to combine the first $N$ linear equations and the last one (i.e. the Kutta-condition) in the Q7 and obtain the matrix form of the whole system of linear equations.

# ----------------------------

# The equations can be solved now! This is the vortex-source panel method.

# --------------------
Please ignore the cell below. It just loads our style for the notebook.
# In[1]:


from IPython.core.display import HTML
def css_styling(filepath):
    styles = open(filepath, 'r').read()
    return HTML(styles)
css_styling('../styles/custom.css')