<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 Content under Creative Commons Attribution license CC-BY 4.0, code under BSD 3-Clause License © 2018 parts of this notebook are from Derivative Approximation by Finite Differences by David Eberly,  additional text and SymPy examples by D. Koehn, notebook style sheet by L.A. Barba, N.C. Clementi
Step1:  Generalization of Taylor FD operators
In the last lesson, we learned how to derive a high order FD approximation for the second derivative using Taylor series expansion. In the next step we derive a general equation to compute FD operators, where I use a detailed derivation based on "Derivative Approximation by Finite Differences" by David Eberly
Estimation of arbitrary FD operators by Taylor series expansion
We can approximate the $d-th$ order derivative of a function $f(x)$ with an order of error $p>0$ by a general finite-difference approximation
<ASSISTANT_TASK:>
Python Code:
# Execute this cell to load the notebook's style sheet, then ignore it
from IPython.core.display import HTML
css_file = '../../style/custom.css'
HTML(open(css_file, "r").read())
Explanation: Content under Creative Commons Attribution license CC-BY 4.0, code under BSD 3-Clause License © 2018 parts of this notebook are from Derivative Approximation by Finite Differences by David Eberly,  additional text and SymPy examples by D. Koehn, notebook style sheet by L.A. Barba, N.C. Clementi
End of explanation
# import SymPy libraries
from sympy import symbols, differentiate_finite, Function
# Define symbols
x, h = symbols('x h')
f = Function('f')
Explanation: Generalization of Taylor FD operators
In the last lesson, we learned how to derive a high order FD approximation for the second derivative using Taylor series expansion. In the next step we derive a general equation to compute FD operators, where I use a detailed derivation based on "Derivative Approximation by Finite Differences" by David Eberly
Estimation of arbitrary FD operators by Taylor series expansion
We can approximate the $d-th$ order derivative of a function $f(x)$ with an order of error $p>0$ by a general finite-difference approximation:
\begin{equation}
\frac{h^d}{d!}f^{(d)}(x) = \sum_{i=i_{min}}^{i_{max}} C_i f(x+ih) + \cal{O}(h^{d+p})
\end{equation}
where h is an equidistant grid point distance. By choosing the extreme indices $i_{min}$ and $i_{max}$, you can define forward, backward or central operators. The accuracy of the FD operator is defined by it's length and therefore also the number of 
weighting coefficients $C_i$ incorporated in the approximation. $\cal{O}(h^{d+p})$ terms are negelected. 
Formally, we can approximate $f(x+ih)$ by a Taylor series expansion:
\begin{equation}
f(x+ih) = \sum_{n=0}^{\infty} i^n \frac{h^n}{n!}f^{(n)}(x)\nonumber
\end{equation}
Inserting into eq.(1) yields
\begin{align}
\frac{h^d}{d!}f^{(d)}(x) &= \sum_{i=i_{min}}^{i_{max}} C_i \sum_{n=0}^{\infty} i^n \frac{h^n}{n!}f^{(n)}(x) + \cal{O}(h^{d+p})\nonumber\
\end{align}
We can move the second sum on the RHS to the front
\begin{align}
\frac{h^d}{d!}f^{(d)}(x) &= \sum_{n=0}^{\infty} \left(\sum_{i=i_{min}}^{i_{max}} i^n C_i\right) \frac{h^n}{n!}f^{(n)}(x) + \cal{O}(h^{d+p})\nonumber\
\end{align}
In the FD approximation we only expand the Taylor series up to the term $n=(d+p)-1$
\begin{align}
\frac{h^d}{d!}f^{(d)}(x) &= \sum_{n=0}^{(d+p)-1} \left(\sum_{i=i_{min}}^{i_{max}} i^n C_i\right) \frac{h^n}{n!}f^{(n)}(x) + \cal{O}(h^{d+p})\nonumber\
\end{align}
and neglect the $\cal{O}(h^{d+p})$ terms
\begin{align}
\frac{h^d}{d!}f^{(d)}(x) &= \sum_{n=0}^{(d+p)-1} \left(\sum_{i=i_{min}}^{i_{max}} i^n C_i\right) \frac{h^n}{n!}f^{(n)}(x)\
\end{align}
Multiplying by $\frac{d!}{h^d}$ leads to the desired approximation for the $d-th$ derivative of the function f(x):
\begin{align}
f^{(d)}(x) &= \frac{d!}{h^d}\sum_{n=0}^{(d+p)-1} \left(\sum_{i=i_{min}}^{i_{max}} i^n C_i\right) \frac{h^n}{n!}f^{(n)}(x)\
\end{align}
Treating the approximation in eq.(2) as an equality, the only term in the sum on the right-hand side of the approximation that contains $\frac{h^d}{d!}f^{d}(x)$ occurs when $n = d$, so the coefficient of that term must be 1. The other terms must vanish for there to be equality, so the coefficients of those terms must be 0; therefore, it is necessary that
\begin{equation}
\sum_{i=i_{min}}^{i_{max}} i^n C_i=
\begin{cases}
0, ~~ 0 \le n \le (d+p)-1 ~ \text{and} ~ n \ne d\
1, ~~ n = d
\end{cases}\nonumber\
\end{equation}
This is a set of $d + p$ linear equations in $i_{max} − i_{min} + 1$ unknowns. If we constrain the number of unknowns to be $d+p$, the linear system has a unique solution. 
A forward difference approximation occurs if we set $i_{min} = 0$
and $i_{max} = d + p − 1$. 
A backward difference approximation can be implemented by setting $i_{max} = 0$ and $i_{min} = −(d + p − 1)$.
A centered difference approximation occurs if we set $i_{max} = −i_{min} = (d + p − 1)/2$ where it appears that $d + p$ is necessarily an odd number. As it turns out, $p$ can be chosen to be even regardless of the parity of $d$ and $i_{max} = (d + p − 1)/2$.
We could either implement the resulting linear system as matrix equation as in the previous lesson, or simply use a SymPy function which gives us the FD operators right away.
End of explanation
<END_TASK>