<jupyter_start><jupyter_text>Computational Physics Project 1 -- Interactive Fourier Serise Demo <NAME> August 27th 2018 For use with python 3.6 or above<jupyter_code>import numpy as np
from scipy.signal import square
import matplotlib.pyplot as plt
from IPython.display import HTML, Image
from ipywidgets import widgets
from scipy.signal import square, sawtooth

plt.rc('text', usetex=True)
plt.rc('font', family='Serif')

%matplotlib inline<jupyter_output><empty_output><jupyter_text>IndroductionThe fourier serise is the most obvious, and analytic, manifestation of Fourier's Theorum (Taylor 192). The serise is defined in its purest form as the sumation of an infinite set of sin waves of varying amplitides and angular frequencies. The general form of the serise if given as:$$f(t) = \sum^{\infty}_{n=0}[a_{n}\cos(n\omega t) + b_{n}\sin(n\omega t)]$$We can expand the serise to $n$ terms to fit to any periodic function to an arbitrary level of presicion.Interactive Fourier Serise Below I setup a set of peicwise defined functions which will be used to demonstrate how a Fourier Serise can be used. I define the following functions:   * $f(x) = \begin{cases}\sin(x) & x  0.5\end{cases}$      * $f(x) = \begin{cases}2x+1 & x < 0 \\ 2x & x \geq 0\end{cases}$along with a standard square wave, and a standard sawtooth wave.<jupyter_code>def func(x, func='square', wavelength=1):
    if func == 'square':
        sig = np.sin(2 * 1/wavelength * np.pi * x)
        return square(2 * 1/wavelength * np.pi * x, duty=(sig + 1)/2)
    if func == 'crazy':
        if x < -0.5:
            return np.sin(x)
        if -0.5 <= x <= 0.5:
            return 0.5
        if x > 0.5:
            return np.sin(x)
    if func == 'saw':
        return sawtooth(2 * 1/wavelength * np.pi * x)
    if func == 'saw2':
        if x < 0:
            return 2*x+1
        else:
            return 2*x<jupyter_output><empty_output><jupyter_text>The coefficients of the Fourier Serise can be calculated as followsFirst Recall that the fourier serise is defined as$$   f(t) = \sum^{\infty}_{n=0}[a_{n}\cos(n\omega t) + b_{n}\sin(n\omega t)]$$We will for the moment assume that all functions we work with are even, such that their behavior can be captured by a $\cos(x(t))$ function. This allows us to focus on finding the values of $a_{n}$. In the case the function is even the fourier serise is defined as$$    f(t) = \sum^{\infty}_{n=0}a_{n}\cos(n\omega t)$$We can multiply by one in the form of $\frac{\cos(x(t))}{\cos(x(t))}\frac{\text{d}t}{\text{d}t}$$$    f(t) = \sum^{\infty}_{n=0}a_{n}\cos(n\omega t)\frac{\cos(m\omega t)}{\cos(m\omega t)}\frac{\text{d}t}{\text{d}t}$$where $m$ is some seperate incriment variable from $n$. Rearanging $$    \cos(m\omega t) f(t)\text{d}t = \sum^{\infty}_{n=0}a_{n}\cos(n\omega t)\cos(m\omega t)\text{d}t$$if we then integrate this over some period$$\int_{T_{0}}^{T}\cos(m\omega t) f(t)\text{d}t = \int_{T_{0}}^{T}\sum^{\infty}_{n=0}a_{n}\cos(n\omega t)\cos(m\omega t)\text{d}t \\\int_{T_{0}}^{T}\cos(m\omega t) f(t)\text{d}t = \sum^{\infty}_{n=0}a_{n}\int_{T_{0}}^{T}\cos(n\omega t)\cos(m\omega t)\text{d}t$$Invoking trigonometric identities we can say $$\cos(n\omega t)\cos(m\omega t) = \frac{\cos((m+n)\omega t)+\cos((m-n)\omega t)}{2}$$so then$$2\int_{T_{0}}^{T}\cos(m\omega t) f(t)\text{d}t = \sum^{\infty}_{n=0}a_{n}\int_{T_{0}}^{T}\cos((m+n)\omega t)+\cos((m-n)\omega t)\text{d}t$$if we define $m=\{m_{i} \in \mathbb{Z}$ | $m_{i}\geq 0\}$ then the integral of the first term in the integrand goes to 0 (as all for all m and n where $n=\{n_{i} \in \mathbb{Z}$ | $-\infty \leq n_{i}\leq \infty\}$ that sine term will cover an integer number of wavelengths and so contain net 0 area). So then we can say:$$2\int_{T_{0}}^{T}\cos(m\omega t) f(t)\text{d}t = \sum^{\infty}_{n=0}a_{n}\int_{T_{0}}^{T}\cos((m-n)\omega t)\text{d}t$$because cosine is an even function$$\int_{T_{0}}^{T}\cos((m-n)\omega t)\text{d}t = 0(1-\delta_{mn}) + (T-T_{0})\delta(mn) \\\int_{T_{0}}^{T}\cos((m-n)\omega t)\text{d}t = (T-T_{0})\delta(mn)$$so then$$2\int_{T_{0}}^{T}\cos(m\omega t) f(t)\text{d}t =  (T-T_{0})\delta(mn)\sum^{\infty}_{n=0}a_{n}$$if we define $T_{0} = 0$ then$$2\int_{T_{0}}^{T}\cos(m\omega t) f(t)\text{d}t =  [T\delta(mn)]\sum^{\infty}_{n=0}a_{n}$$but because m must equal n for a non--trivial solution$$2\int_{T_{0}}^{T}\cos(n\omega t) f(t)\text{d}t =  Ta_{n}$$which can be rearanged too$$a_{n} = \frac{2}{T}\int_{T_{0}}^{T}\cos(n\omega t) f(t)\text{d}t$$A similar set of steps can be followed to derive the coefficients $b_{n}$. This full derivation can be found here. Below I define two functions which calculate these coefficients (using a quire naive numerical integration method). I have refactored these in terms of wavelength and wavenumber (as opposed to period and angular frequency), this is because I think better in wavelength space than in period space.<jupyter_code>def A_m(m, xvalues, dx, k, wavelength, ftype='square'):
    s = 0
    for x in xvalues:
        s += func(x, func=ftype)*np.cos(m*k*x)*dx
    return (2/wavelength)*s
def B_m(m, xvalues, dx, k, wavelength, ftype='square'):
    s = 0
    for x in xvalues:
        s += func(x, func=ftype)*np.sin(m*k*x)*dx
    return (2/wavelength)*s<jupyter_output><empty_output><jupyter_text>I now define a function which will be used to interactivly plot any of the functions defined in the function > func(x, func)<jupyter_code>def main(N=4000, wavelength=1, numTerms=6, ftype='square', returnFig=False):
    k=(2*np.pi/wavelength)
    dx=wavelength/N
    integers = np.arange(1, numTerms+1)
    xvalues = np.linspace(-wavelength/2, wavelength/2, N)
    xplot = np.linspace(-wavelength, wavelength, N)
    
    A_terms = np.zeros(N)
    B_terms = np.zeros(N)
    f_of_x = np.zeros(N)
    
    A0_term = 0.5*A_m(0, xvalues, dx, k, wavelength, ftype=ftype)
    
    for m in integers:
        A_terms += A_m(m, xvalues, dx, k, wavelength, ftype=ftype)*np.cos(m*k*xplot)
        B_terms += B_m(m, xvalues, dx, k, wavelength, ftype=ftype)*np.sin(m*k*xplot)
        
    f_of_x = A0_term + A_terms + B_terms
    true = [func(i, func=ftype) for i in xplot]
    
    fig = plt.figure(figsize=(10, 7))
    signal = fig.add_subplot(211)
    signal.plot(xplot, true, label='True')
    signal.plot(xplot, f_of_x, label='predicted')
    signal.set_ylabel('Dependant', fontsize=20)
    signal.set_xlabel('Independat', fontsize=20)
    signal.tick_params(which='both', labelsize=15)
    plt.legend(loc='best')
    
    residuals = fig.add_subplot(212)
    residuals.plot(xplot, [x-y for x, y in zip(true, f_of_x)], 'o')
    residuals.set_ylabel('Data-Serise', fontsize=20)
    residuals.set_xlabel('Independat', fontsize=20)
    
    residuals.tick_params(which='both',labelsize=15)

    if returnFig:
        return fig, (signal, residuals)
widgets.interact(main, numTerms=(1, 100), N=(10, 500), wavelength=widgets.fixed(2), returnFig=widgets.fixed(False))<jupyter_output><empty_output><jupyter_text>One can adjust the number of terms in the serise, as well as the number of points used in the plotting of the serise. One can also enter the name of the function (as defined in the function func(x, func)) to change which function the fourier serise is plotting.By increasing the number of terms we see that the serise better approximated the data (however there will always be some strange action at any disconituites in the fit data). Conversly we see that fewer terms fit the data poorly. One interesting to note is the diminishing returnes in terms of accuracy of serise as number of terms terms increase. That is to say that going from 10 $\rightarrow$ 20 terms makes a much alrger difference (visually at least) than going from say 500 $\rightarrow$ 600 terms does.Also note which functions the fourier serise manages to fit. Periodic functions such as my peicwise defined square wave, or a sawtooth wave are good candidates for fitting with a fourier serise. However, while non periodic functions can still be fit, the serise generally needs more terms to capture their structure, and obviously cannot make predictive statements about their long term behavior if they don't exibit any periodicities. Animation We can use the python module mplEasyAnimate to make a simple video which we can embed. Below we make a video demonstrating a steady increase in the number of terms of the fourier Serise.<jupyter_code>numTerms = np.linspace(1, 100, 100)
N = 500
from mplEasyAnimate import animation
from tqdm import tqdm, tqdm_notebook

anim = animation('FTermIncrease.mp4')

for numTerm in tqdm_notebook(numTerms):
    fig, ax = main(N=4000, wavelength=1, numTerms=numTerm, ftype='square', returnFig=True)
    ax[0].annotate(f'Terms: {numTerm}', xy=(-0.5, 0.5), fontsize=17)
    anim.add_frame(fig)
    plt.close(fig)
anim.close()
import io
import base64
from IPython.display import HTML
import os

if not os.path.exists('FTermIncrease.mp4'):
    raise IOError('ERROR! Animation has not been generated to the local directory yet!')
video = io.open('FTermIncrease.mp4', 'r+b').read()
encoded = base64.b64encode(video)
HTML(data='''<video alt="test" controls>
                <source src="data:video/mp4;base64,{0}" type="video/mp4" />
             </video>'''.format(encoded.decode('ascii')))<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1491331969419619, 'fraction_numerical': 0.019597286529557445, 'mean_word_length': 4.177257525083612, 'pattern_counts': {'":': 0, '<': 30, '<?xml version=': 0, '>': 27, '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': '5', 'dir': 'jupyter-structured-clean-dedup', 'id': '178171', 'license': 'MIT', 'n_tokens_mistral': 3370, 'n_tokens_neox': 3208, 'n_words': 1005, 'path': 'InteractiveFourierSerise/BoudreauxProject.ipynb', 'repo_name': 'tboudreaux/ComputationalPhysics', 'size': '9300'}