<jupyter_start><jupyter_text>Finite Difference MethodThis note book illustrates the finite different method for a Boundary Value Problem. Example Boudary Value Problem$$  \frac{d^2 y}{dx^2} = 4y$$ Boundary Condition$$ y(0)=1.1752,  y(1)=10.0179 $$<jupyter_code>import numpy as np
import math
import matplotlib.pyplot as plt
import warnings
warnings.filterwarnings("ignore")
class ListTable(list):
    """ Overridden list class which takes a 2-dimensional list of 
        the form [[1,2,3],[4,5,6]], and renders an HTML Table in 
        IPython Notebook. """
from IPython.core.display import HTML

def _repr_html_(self):
    html = ["<table>"]
    for row in self:
        html.append("<tr>")

        for col in row:
            html.append("<td>{0}</td>".format(col))

        html.append("</tr>")
    html.append("</table>")
    return ''.join(html)
from IPython.display import HTML    
HTML('''<script>
code_show=true; 
function code_toggle() {
 if (code_show){
 $('div.input').hide();
 } else {
 $('div.input').show();
 }
 code_show = !code_show
} 
$( document ).ready(code_toggle);
</script>
<form action="javascript:code_toggle()"><input type="submit" value="Click here to toggle on/off the raw code."></form>''')<jupyter_output><empty_output><jupyter_text>Discrete AxisThe stepsize is defined as$$h=\frac{b-a}{N}$$here it is $$h=\frac{1-0}{10}$$giving $$x_i=0+0.1 i$$for $i=0,1,...10.$<jupyter_code>## BVP
N=10
h=1/N
x=np.linspace(0,1,N+1)
fig = plt.figure(figsize=(10,4))
plt.plot(x,0*x,'o:',color='red')
plt.xlim((0,1))
plt.xlabel('x',fontsize=16)
plt.title('Illustration of discrete time points for h=%s'%(h),fontsize=32)
plt.show()<jupyter_output><empty_output><jupyter_text>The Difference EquationThe gerenal difference equation is$$ \frac{1}{h^2}\left(y_{i-1}-2y_i+y_{i+1}\right)=4y_i \ \ \ i=1,..,N-1. $$Rearranging the equation we have the system of N-1 equations$$i=1: \frac{1}{0.1^2}\color{green}{y_{0}} -\left(\frac{2}{0.1^2}+4\right)y_1 +\frac{1}{0.1^2} y_{2}=0$$$$i=2: \frac{1}{0.1^2}y_{1} -\left(\frac{2}{0.1^2}+4\right)y_2 +\frac{1}{0.1^2} y_{3}=0$$$$ ...$$$$i=8: \frac{1}{0.1^2}y_{7} -\left(\frac{2}{0.1^2}+4\right)y_8 +\frac{1}{0.1^2} y_{9}=0$$$$i=9: \frac{1}{0.1^2}y_{8} -\left(\frac{2}{0.1^2}+4\right)y_9 +\frac{1}{0.1^2} \color{green}{y_{10}}=0$$where the green terms are the known boundary conditions.Rearranging the equation we have the system of 9 equations$$i=1:  -\left(\frac{2}{0.1^2}+4\right)y_1 +\frac{1}{0.1^2} y_{2}=0$$$$i=2: \frac{1}{0.1^2}y_{1} -\left(\frac{2}{0.1^2}+4\right)y_2 +\frac{1}{0.1^2} y_{3}=-\frac{1}{0.1^2}\color{green}{y_{0}}$$$$ ...$$$$i=8: \frac{1}{0.1^2}y_{7} -\left(\frac{2}{0.1^2}+4\right)y_8 +\frac{1}{0.1^2} y_{9}=0$$$$i=9: \frac{1}{0.1^2}y_{8} -\left(\frac{2}{0.1^2}+4\right)y_9 =0$$where the green terms are the known boundary conditions.Putting this into matrix form gives a $9\times 9 $ matrix $$A=\left(\begin{array}{ccc ccc ccc} -204&100&0& 0&0&0& 0&0&0\\100&-204&100 &0&0&0& 0&0&0\\0&100&-204& 100&0&0& 0&0&0\\.&.&.& .&.&.& .&.&.\\.&.&.& .&.&.& .&.&.\\0&0&0& 0&0&0& 100&-204&100\\0&0&0& 0&0&0& 0&100&-204\end{array}\right)$$an unknown vector $$\color{red}{\mathbf{y}}=\color{red}{\left(\begin{array}{c} y_1\\y_2\\y_3\\.\\.\\y_8\\y_9\end{array}\right)}$$<jupyter_code>y=np.zeros((N+1))
# Boundary Condition
y[0]=1.1752
y[N]=10.0179<jupyter_output><empty_output><jupyter_text>and the known right hand side is a known $9\times 1$ vector with the boundary conditions$$\mathbf{b}=\left(\begin{array}{c}-117.52\\0\\0\\.\\.\\0\\-1001.79 \end{array}\right)$$$$ A\mathbf{y}=\mathbf{b}$$The plot below is a graphical representation of the matrix A.<jupyter_code>b=np.zeros(N-1)

# Boundary Condition
b[0]=-y[0]/(h*h)
b[N-2]=-y[N]/(h*h)

A=np.zeros((N-1,N-1))
# Diagonal
for i in range (0,N-1):
    A[i,i]=-(2/(h*h)+4)

for i in range (0,N-2):           
    A[i+1,i]=1/(h*h)
    A[i,i+1]=1/(h*h)
    
plt.imshow(A)
plt.xlabel('i',fontsize=16)
plt.ylabel('j',fontsize=16)
plt.yticks(np.arange(N-1), np.arange(1,N-0.9,1))
plt.xticks(np.arange(N-1), np.arange(1,N-0.9,1))
clb=plt.colorbar()
clb.set_label('Matrix value')
plt.title('Matrix A',fontsize=32)
plt.tight_layout()
plt.subplots_adjust()
plt.show()<jupyter_output><empty_output><jupyter_text>Solving the systemTo solve invert the matrix $A$ such that $$A^{-1}Ay=A^{-1}b$$$$y=A^{-1}b$$The plot below shows the graphical representation of $A^{-1}$.<jupyter_code>invA=np.linalg.inv(A)

plt.imshow(invA)
plt.xlabel('i',fontsize=16)
plt.ylabel('j',fontsize=16)
plt.yticks(np.arange(N-1), np.arange(1,N-0.9,1))
plt.xticks(np.arange(N-1), np.arange(1,N-0.9,1))
clb=plt.colorbar()
clb.set_label('Matrix value')
plt.title(r'Matrix $A^{-1}$',fontsize=32)
plt.tight_layout()
plt.subplots_adjust()
plt.show()


y[1:N]=np.dot(invA,b)<jupyter_output><empty_output><jupyter_text>Result The plot below shows the approximate solution of the Boundary Value Problem (blue v) and the exact solution (black dashed line).<jupyter_code>fig = plt.figure(figsize=(8,4))

plt.plot(x,y,'v',label='Finite Difference')
plt.plot(x,np.sinh(2*x+1),'k:',label='exact')
plt.xlabel('x')
plt.ylabel('y')
plt.legend(loc='best')
plt.show()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.24489795918367346, 'fraction_numerical': 0.07094266277939747, 'mean_word_length': 5.207478890229192, 'pattern_counts': {'":': 0, '<': 36, '<?xml version=': 0, '>': 36, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '6', 'dir': 'jupyter-structured-clean-dedup', 'id': '13278', 'license': 'MIT', 'n_tokens_mistral': 2400, 'n_tokens_neox': 2242, 'n_words': 371, 'path': 'Chapter 06 - Boundary Value Problems/.ipynb_checkpoints/Boundary Value Problem-checkpoint.ipynb', 'repo_name': 'jjcrofts77/Numerical-Analysis-Python', 'size': '5145'}