
!pip3 install numpy scipy matplotlib uncertainties

Defaulting to user installation because normal site-packages is not writeable
Requirement already satisfied: numpy in /home/nicole/.local/lib/python3.9/site-packages (1.19.5)
Requirement already satisfied: scipy in /home/nicole/.local/lib/python3.9/site-packages (1.6.0)
Requirement already satisfied: matplotlib in /home/nicole/.local/lib/python3.9/site-packages (3.3.4)
Requirement already satisfied: uncertainties in /home/nicole/.local/lib/python3.9/site-packages (3.1.5)
Requirement already satisfied: cycler>=0.10 in /home/nicole/.local/lib/python3.9/site-packages (from matplotlib) (0.10.0)
Requirement already satisfied: numpy in /home/nicole/.local/lib/python3.9/site-packages (1.19.5)
Requirement already satisfied: pyparsing!=2.0.4,!=2.1.2,!=2.1.6,>=2.0.3 in /usr/lib/python3.9/site-packages (from matplotlib) (2.4.7)
Requirement already satisfied: pillow>=6.2.0 in /usr/lib/python3.9/site-packages (from matplotlib) (8.1.0)
Requirement already satisfied: kiwisolver>=1.0.1 in /home/nicole/.local/lib/python3.9/site-packages (from matplotlib) (1.3.1)
Requirement already satisfied: python-dateutil>=2.1 in /usr/lib/python3.9/site-packages (from matplotlib) (2.8.1)
Requirement already satisfied: numpy in /home/nicole/.local/lib/python3.9/site-packages (1.19.5)
Requirement already satisfied: future in /home/nicole/.local/lib/python3.9/site-packages (from uncertainties) (0.18.2)
Requirement already satisfied: six in /usr/lib/python3.9/site-packages (from cycler>=0.10->matplotlib) (1.15.0)
Requirement already satisfied: six in /usr/lib/python3.9/site-packages (from cycler>=0.10->matplotlib) (1.15.0)


# from google.colab import drive
# drive.mount('/content/drive')


# %cd '/content/drive/MyDrive/Seminaire python 25-05-2021'
!ls  # check you're in the correct directory

 assets
 data
'demo0 basics.py'
'demo1 fib print.py'
'demo2 fib list.py'
'demo3 (advanced) fib yield.py'
'demo4 fib open.py'
'demo5 (advanced) fib open list comprehension.py'
'demo6 numpy matplotlib matlab-style.py'
'demo7 numpy matplotlib python-style.py'
'demo8 notebook.html'
'demo8 notebook.ipynb'
 generated
 presentation.pdf
 presentation.pptx
 resources


# !apt install texlive texlive-latex-extra texlive-fonts-recommended dvipng cm-super


import numpy as np
import matplotlib.pyplot as plt
import scipy.constants as cst


a = 1.05  # arb. units
lamb0 = 643.8*1e-9  # m
d = 4.04*1e-3  # m
n = 1.4567  # no units


data_mB, data_2s = np.loadtxt('data/interferometer/parallel.csv', delimiter=',', skiprows=1, unpack=True)


data_B = 1e-3*data_mB  # convert to T
data_s = data_2s/2
rhs = cst.h*cst.c * data_s/a * 1/(2*d) * 1/np.sqrt(n**2-1)

# fig, ax = plt.subplots()
fig, ax = plt.subplots(dpi=120, facecolor='white')
ax.plot(data_B, rhs, 'o')
ax.set_xlabel(r'$B$ [T]')
ax.set_ylabel(r'$\frac{s}{a} \frac{hc}{2d} \frac{1}{\sqrt{n^2-1}}$ [J]')
plt.show(fig)
plt.close(fig)


fit_coefs = np.polyfit(data_B, rhs, 1)
print(fit_coefs)

[8.67602362e-24 8.67302103e-25]


fit_coefs, fit_cov = np.polyfit(data_B, rhs, 1, cov=True)
print(fit_cov)

[[ 3.96666320e-50 -2.04917821e-50]
 [-2.04917821e-50  1.10927102e-50]]


fit_stddevs = np.sqrt(np.diag(fit_cov))
print(fit_stddevs)

[1.99164836e-25 1.05321936e-25]


mu_B_exp = fit_coefs[0]
mu_B_exp_err = fit_stddevs[0]
print('mu_B = ' + str(mu_B_exp) + ' +/- ' + str(mu_B_exp_err))  # matlab-like
print('mu_B = %.3e +/- %.3e' % (mu_B_exp, mu_B_exp_err))  # c printf-like
print('mu_B = {:.3e} +/- {:.3e}'.format(mu_B_exp, mu_B_exp_err))  # .format method
print(f'mu_B = {mu_B_exp:.3e} +/- {mu_B_exp_err:.3e}')  # format string

mu_B = 8.676023621077322e-24 +/- 1.9916483634199414e-25
mu_B = 8.676e-24 +/- 1.992e-25
mu_B = 8.676e-24 +/- 1.992e-25
mu_B = 8.676e-24 +/- 1.992e-25


print('{:.1%}'.format(abs(mu_B_exp-cst.value('Bohr magneton'))/cst.value('Bohr magneton')))

6.4%


def SS_res(x, y):
    fit_coefs = np.polyfit(x, y, 1)
    f = np.poly1d(fit_coefs)
    return np.sum((y - f(x))**2)

def SS_tot(y):
    return np.sum((y - np.mean(y))**2)

def R2(x, y):
    return 1 - SS_res(x, y)/SS_tot(y)

print(SS_res(data_B, rhs))  # our SS_res
print(np.polyfit(data_B, rhs, 1, full=True)[1])  # numpy.polyfit residual

print('our R² ≈ {:.5f}'.format(R2(data_B, rhs)))
import scipy.stats
print('scipy.stats R² ≈ {:.5f}'.format(scipy.stats.linregress(data_B, rhs)[2]**2))

4.0532443545891705e-50
[4.05324435e-50]
our R² ≈ 0.99580
scipy.stats R² ≈ 0.99580


def Pearson_r(x, y):
    a = np.sum((x - np.mean(x))*(y - np.mean(y)))
    b = np.sqrt(np.sum((x - np.mean(x))**2)) * np.sqrt(np.sum((y - np.mean(y))**2))
    return a/b

print('our r ≈ {:.5f}'.format(Pearson_r(data_B, rhs)))
print('scipy.stats.pearsonr r ≈ {:.5f}'.format(scipy.stats.pearsonr(data_B, rhs)[0]))

our r ≈ 0.99790
scipy.stats.pearsonr r ≈ 0.99790


# r² = R2 up to machine precision
print(Pearson_r(data_B, rhs)**2 - R2(data_B, rhs))

5.551115123125783e-16


fit_poly = np.poly1d(fit_coefs)

fig, ax = plt.subplots(dpi=120, facecolor='white')
ax.plot(data_B, rhs, 'o')
# ax.plot(data_B, fit_poly(data_B), '--')
ax.plot(
    data_B, fit_poly(data_B),
    '--', color='tab:grey',
    label=f'Fit for $B_\parallel$ $\\sim {fit_coefs[0]:.2e} \pm {fit_stddevs[0]:.2e}B$'
)
ax.legend()
ax.set_xlabel(r'$B$ [T]')
ax.set_ylabel(r'$\frac{s}{a} \frac{hc}{2d} \frac{1}{\sqrt{n^2-1}}$ [J]')
plt.show(fig)
plt.close(fig)


from uncertainties import ufloat  # container for nominal value + error
import uncertainties.unumpy as unp  # math with uncertainty propagation


u_a = ufloat(1.050, 0.025)  # arb units
u_lamb0 = ufloat(643.8, 0.05)*1e-9  # m
u_d = ufloat(4.04, 0.005)*1e-3  # m
u_n = ufloat(1.4567, 0.00005)  # no units
print(u_a, u_lamb0, u_d, u_n, sep='\n')

1.050+/-0.025
(6.4380+/-0.0005)e-07
0.004040+/-0.000005
1.45670+/-0.00005


u_B = unp.uarray(data_B, 1e-3)
u_s = unp.uarray(data_s, 0.025)

u_rhs = cst.h*cst.c * u_s/u_a * 1/(2*u_d) * 1/unp.sqrt(u_n**2-1)


# This doesn't work !
# fit_coefs, fit_cov = np.polyfit(u_B, u_rhs, 1, cov=True)


fit_coefs, fit_cov = np.polyfit(
    unp.nominal_values(u_B), unp.nominal_values(u_rhs),
    1, cov=True
)
fit_stddevs = np.sqrt(np.diag(fit_cov))
print(fit_coefs, fit_stddevs, sep='\n')

[8.67602362e-24 8.67302103e-25]
[1.99164836e-25 1.05321936e-25]


u_fit_coefs = unp.uarray(fit_coefs, fit_stddevs)
u_mu_B = u_fit_coefs[0]
# or
# u_mu_B = ufloat(fit_coefs[0], fit_stddevs[0])
print(u_mu_B)
print(f'{u_mu_B:.2L}')  # LaTeX formatting

(8.68+/-0.20)e-24
\left(8.7 \pm 0.2\right) \times 10^{-24}


fit_poly = np.poly1d(fit_coefs)

fig, ax = plt.subplots(dpi=120, facecolor='white')
ax.plot(
    unp.nominal_values(u_B), fit_poly(unp.nominal_values(u_B)),
    '--', color='tab:grey',
    label=f'Fit for $B_\\parallel$ $\\sim {u_fit_coefs[0]:.2L}B$'
)
ax.errorbar(
    unp.nominal_values(u_B), unp.nominal_values(u_rhs),
    xerr=unp.std_devs(u_B), yerr=unp.std_devs(u_rhs),
    marker='.', linestyle=''
)
ax.legend(loc='upper left')
ax.set_xlabel(r'$B$ [T]')
ax.set_ylabel(r'$\frac{s}{a} \frac{hc}{2d} \frac{1}{\sqrt{n^2-1}}$ [J]')
plt.show(fig)
plt.close(fig)


import glob
import os.path

print(glob.glob('data/interferometer/*.csv'))
print([os.path.basename(filepath) for filepath in glob.glob('data/interferometer/*.csv')])

['data/interferometer/transverse.csv', 'data/interferometer/parallel.csv']
['transverse.csv', 'parallel.csv']


fig, ax = plt.subplots(dpi=120, facecolor='white')

for filepath in glob.glob('data/interferometer/*.csv'):
    # Load the data
    data_B, data_2s = np.loadtxt(filepath, delimiter=',', skiprows=1, unpack=True)
    u_B = unp.uarray(data_B, 2)*1e-3  # T
    u_s = unp.uarray(data_2s, 0.05)/2  # arb units
    
    # Compute the rhs
    u_rhs = cst.h*cst.c * u_s/u_a * 1/(2*u_d) * 1/unp.sqrt(u_n**2-1)
    
    # Perform the fit
    fit_coefs, fit_cov = np.polyfit(unp.nominal_values(u_B), unp.nominal_values(u_rhs), 1, cov=True)
    fit_ucoefs = unp.uarray(fit_coefs, np.sqrt(np.diag(fit_cov)))
    fit_poly = np.poly1d(fit_coefs)
    
    # Plot stuff
    line, = ax.plot(
        unp.nominal_values(u_B), fit_poly(unp.nominal_values(u_B)),
        '--',
        label=f'Fit for $B$ $\\sim {fit_ucoefs[0]:L}B$'
    )
    ax.errorbar(
        unp.nominal_values(u_B), unp.nominal_values(u_rhs),
        xerr=unp.std_devs(u_B), yerr=unp.std_devs(u_rhs),
        marker='.', linestyle='', color=line.get_color(),
    )

ax.set_xlabel(r'$B$ [T]')
ax.set_ylabel(r'$\frac{s}{a} \frac{hc}{2d} \frac{1}{\sqrt{n^2-1}}$ [J]')
ax.legend()
fig.savefig('generated/interferometer.pdf')
plt.show(fig)
plt.close(fig)


symbol_lookup = {
    'transverse.csv': '\\perp',
    'parallel.csv': '\\parallel'
}
print(symbol_lookup.get('transverse.csv'))
print(symbol_lookup.get('asdf.csv'))
print(symbol_lookup.get('asdf.csv', 'FALLBACK'))
# or
# print(symbol_lookup['transverse.csv'])
# works too, but doesn't have a fallback option

\perp
None
FALLBACK


fig, ax = plt.subplots(dpi=120, facecolor='white')

for filepath in glob.glob('data/interferometer/*.csv'):
    # Load the data
    data_B, data_2s = np.loadtxt(filepath, delimiter=',', skiprows=1, unpack=True)
    u_B = unp.uarray(data_B, 2)*1e-3  # T
    u_s = unp.uarray(data_2s, 0.05)/2  # arb units
    
    # Compute the rhs
    u_rhs = cst.h*cst.c * u_s/u_a * 1/(2*u_d) * 1/unp.sqrt(u_n**2-1)
    
    # Perform the fit
    fit_coefs, fit_cov = np.polyfit(unp.nominal_values(u_B), unp.nominal_values(u_rhs), 1, cov=True)
    fit_ucoefs = unp.uarray(fit_coefs, np.sqrt(np.diag(fit_cov)))
    fit_poly = np.poly1d(fit_coefs)
    
    # Plot stuff
    symbol = symbol_lookup.get(os.path.basename(filepath))
    line, = ax.plot(
        unp.nominal_values(u_B), fit_poly(unp.nominal_values(u_B)),
        '--',
        label=f'Fit for $B_{symbol}$ $\\sim {fit_ucoefs[0]:L}B$'
    )
    ax.errorbar(
        unp.nominal_values(u_B), unp.nominal_values(u_rhs),
        xerr=unp.std_devs(u_B), yerr=unp.std_devs(u_rhs),
        marker='.', linestyle='', color=line.get_color(),
    )

ax.set_xlabel(r'$B$ [T]')
ax.set_ylabel(r'$\frac{s}{a} \frac{hc}{2d} \frac{1}{\sqrt{n^2-1}}$ [J]')
ax.legend()
fig.savefig('generated/interferometer.pdf')
plt.show(fig)
plt.close(fig)


import matplotlib

matplotlib.rcParams['figure.figsize'] = [4.0, 3.0]  # Empirical size for half-page in report
matplotlib.rcParams['font.size'] = 13  # Empirical size for half-page in report
matplotlib.rcParams['text.usetex'] = True  # LaTeX rendering
matplotlib.rcParams['font.family'] = 'serif'  # LaTeX serif font
matplotlib.rcParams['savefig.bbox'] = 'tight'  # Automatically optimize margins


fig, ax = plt.subplots(dpi=150, facecolor='white')

for filepath in glob.glob('data/interferometer/*.csv'):
    # Load the data
    data_B, data_2s = np.loadtxt(filepath, delimiter=',', skiprows=1, unpack=True)
    u_B = unp.uarray(data_B, 2)*1e-3  # T
    u_s = unp.uarray(data_2s, 0.05)/2  # arb units
    
    # Compute the rhs
    u_rhs = cst.h*cst.c * u_s/u_a * 1/(2*u_d) * 1/unp.sqrt(u_n**2-1)
    
    # Perform the fit
    fit_coefs, fit_cov = np.polyfit(unp.nominal_values(u_B), unp.nominal_values(u_rhs), 1, cov=True)
    fit_ucoefs = unp.uarray(fit_coefs, np.sqrt(np.diag(fit_cov)))
    fit_poly = np.poly1d(fit_coefs)
    
    # Plot stuff
    symbol = symbol_lookup.get(os.path.basename(filepath))
    line, = ax.plot(
        unp.nominal_values(u_B), fit_poly(unp.nominal_values(u_B)),
        '--',
        label=f'Fit for $B_{symbol}$ $\\sim {fit_ucoefs[0]:L}B$'
    )
    ax.errorbar(
        unp.nominal_values(u_B), unp.nominal_values(u_rhs),
        xerr=unp.std_devs(u_B), yerr=unp.std_devs(u_rhs),
        marker='.', linestyle='', color=line.get_color(),
    )

ax.set_xlabel(r'$B$ [T]')
ax.set_ylabel(r'$\frac{s}{a} \frac{hc}{2d} \frac{1}{\sqrt{n^2-1}}$ [J]')
ax.legend(fontsize='x-small')
fig.savefig('generated/interferometer.pdf')
plt.show(fig)
plt.close(fig)


data_U, data_I = np.loadtxt(
    'data/plasma/f=70 amp=5.0 p=1.0e-3 d=8.0 I_fil=44 U_grid=0 which=p.csv',
    delimiter=',', unpack=True, skiprows=1
)
data_I *= 1e-4  # Convert to A


fig, ax = plt.subplots(dpi=150, facecolor='white')
ax.plot(data_U, data_I, 'o', markersize=0.5)
ax.set_xlabel('$U$ [V]')
ax.set_ylabel('$I$ [A]')
ax.set_yscale('log')  # optional `base` kwarg, for log2 or ln
plt.show(fig)
plt.close(fig)


from scipy.ndimage.filters import gaussian_filter1d

sigma = 32
U_smooth = gaussian_filter1d(data_U, sigma)
I_smooth = gaussian_filter1d(data_I, sigma)


fig, ax = plt.subplots(dpi=150, facecolor='white')
ax.plot(
    data_U, data_I, 'o',
    markersize=0.5, color='tab:grey', label='Raw data'
)
ax.plot(
    U_smooth, I_smooth,
    color='black', label='Smoothened'
)
ax.set_xlabel('$U$ [V]')
ax.set_ylabel('$I$ [A]')
ax.set_yscale('log')
ax.legend()
plt.show(fig)
plt.close(fig)


arg_from = np.argmin(U_smooth)
arg_to = np.argmax(U_smooth)
U_smooth_trim = U_smooth[arg_from:arg_to]
I_smooth_trim = I_smooth[arg_from:arg_to]


fig, ax = plt.subplots(dpi=150, facecolor='white')
ax.plot(
    data_U, data_I, 'o',
    markersize=0.5, color='tab:grey', label='Raw data'
)
ax.plot(
    U_smooth_trim, I_smooth_trim,
    color='black', label='Smoothened'
)
ax.set_xlabel('$U$ [V]')
ax.set_ylabel('$I$ [A]')
ax.set_yscale('log')
ax.legend()
plt.show(fig)
fig.savefig('generated/plasma.pdf')
plt.close(fig)


import re

# Parse the string
filepath = 'data/plasma/f=70 amp=5.0 p=1.0e-3 d=8.0 I_fil=44 U_grid=0 which=p.csv'
expression = r'(?:f=)(?P<f>[0-9e\.-]+)(?: amp=)(?P<amp>[0-9e\.-]+)(?: p=)(?P<p>[0-9e\.-]+)(?: d=)(?P<d>[0-9e\.-]+)(?: I_fil=)(?P<I_fil>[0-9e\.-]+)(?: U_grid=)(?P<U_grid>[0-9e\.-]+)(?: which=)(?P<which>[pt])(?:\.csv)'
m = re.match(expression, os.path.basename(filepath))
d = m.groupdict()
print(d)

# Convert to float
for key in d.keys():
    if key not in ['which']:
    # or (but not extendable to multiple ignored parameters)
    # if key != 'which'
        d[key] = float(d[key])
print(d)

{'f': '70', 'amp': '5.0', 'p': '1.0e-3', 'd': '8.0', 'I_fil': '44', 'U_grid': '0', 'which': 'p'}
{'f': 70.0, 'amp': 5.0, 'p': 0.001, 'd': 8.0, 'I_fil': 44.0, 'U_grid': 0.0, 'which': 'p'}


x = np.linspace(1, 6, 6)
mask = [True, False, True, True, True, False]
print(x, x[mask], sep='\n')

[1. 2. 3. 4. 5. 6.]
[1. 3. 4. 5.]


print(x % 2 == 0)

[False  True False  True False  True]


print(x[x % 2 == 0])
print(x[x < 4])
print(x[(x % 2 == 0) & (x > 3)])
print(x[(x % 2 == 0) | (x > 3)])

[2. 4. 6.]
[1. 2. 3.]
[4. 6.]
[2. 4. 5. 6.]


fig, ax = plt.subplots(dpi=150, facecolor='white')

for filepath in glob.glob('data/pendulum/initcond *.csv'):
    t, theta, *_ = np.loadtxt(filepath, delimiter=',', unpack=True)
    
    mask = (0 <= t) & (t <= 50)
    t, theta = t[mask], theta[mask]
    
    ax.plot(t, theta)

ax.set_xlabel('$t$ [s]')
ax.set_ylabel('$\\theta$ [arb. units]')
# ax.set_xlim(0, 50)
ax.yaxis.set_major_locator(plt.NullLocator())  # disable yticks
plt.show(fig)
plt.close(fig)


fig, ax = plt.subplots(dpi=150, facecolor='white')

for filepath in glob.glob('data/pendulum/initcond *.csv'):
    t, theta, *_ = np.loadtxt(filepath, delimiter=',', unpack=True)
    
    # Sync up the signals
    dtheta = np.gradient(theta, t)
    arg_start = np.argmax(dtheta < -0.15)
    print(f'{os.path.basename(filepath)} starts at index {arg_start}')
    t, theta = t[arg_start:], theta[arg_start:]
    
    # Apply mask
    mask = (0 <= t) & (t <= 50)
    t, theta = t[mask], theta[mask]
    
    ax.plot(t, theta)

ax.set_xlabel('$t$ [s]')
ax.set_ylabel('$\\theta$ [arb. units]')
ax.yaxis.set_major_locator(plt.NullLocator())  # disable yticks
plt.show(fig)
plt.close(fig)

initcond 300mA 4.9V 2.csv starts at index 38
initcond 300mA 4.9V 4.csv starts at index 42
initcond 300mA 4.9V 1.csv starts at index 28
initcond 300mA 4.9V 3.csv starts at index 30


fig, ax = plt.subplots(dpi=150, facecolor='white')

for filepath in glob.glob('data/pendulum/initcond *.csv'):
    t, theta, *_ = np.loadtxt(filepath, delimiter=',', unpack=True)
    
    # Sync up the signals
    dtheta = np.gradient(theta, t)
    arg_start = np.argmax(dtheta < -0.15)
    t, theta = t[arg_start:], theta[arg_start:]
    
    # Offset time to start at t=0
    t -= t[0]
    
    # Apply mask
    mask = (0 <= t) & (t <= 50)
    t, theta = t[mask], theta[mask]
    
    ax.plot(t, theta)

ax.set_xlabel('$t$ [s]')
ax.set_ylabel('$\\theta$ [arb. units]')
ax.yaxis.set_major_locator(plt.NullLocator())  # disable yticks
plt.show(fig)
fig.savefig('generated/pendulum chaos.pdf')
plt.close(fig)


for filepath in glob.glob('data/pendulum/underdamped *.csv'):
    print(filepath)
    
    t, theta, *_ = np.loadtxt(filepath, delimiter=',', unpack=True)
    
    fig, ax = plt.subplots(ncols=2, dpi=150, facecolor='white', figsize=(8, 3))
    
    ax[0].plot(theta, np.gradient(theta, t), 'o', markersize=0.5)
    ax[0].set_xlabel('$\\theta$')
    ax[0].set_ylabel('$\\dot \\theta$')
    ax[0].xaxis.set_major_locator(plt.NullLocator())
    ax[0].yaxis.set_major_locator(plt.NullLocator())
    
    ax[1].plot(
        np.fft.rfftfreq(len(t), d=np.diff(t)[0]),
        np.abs(np.fft.rfft(theta)),
        linewidth=0.5
    )
    ax[1].set_xlabel('$\\nu$ [Hz]')
    ax[1].set_ylabel('$\mathcal F[\\theta](\\nu)$')
    ax[1].yaxis.set_major_locator(plt.NullLocator())
    ax[1].yaxis.set_minor_locator(plt.NullLocator())
    
    plt.show(fig)
    plt.close(fig)

data/pendulum/underdamped 300mA 4.4V.csv

data/pendulum/underdamped 300mA 4.5V.csv

data/pendulum/underdamped 300mA 5.5V.csv


for filepath in glob.glob('data/pendulum/underdamped *.csv'):
    print(filepath)
    
    t, theta, *_ = np.loadtxt(filepath, delimiter=',', unpack=True)
    
    fig, ax = plt.subplots(ncols=2, dpi=150, facecolor='white', figsize=(8, 3))
    
    ax[0].plot(theta, np.gradient(theta, t), 'o', markersize=0.5)
    ax[0].set_xlabel('$\\theta$')
    ax[0].set_ylabel('$\\dot \\theta$')
    ax[0].xaxis.set_major_locator(plt.NullLocator())
    ax[0].yaxis.set_major_locator(plt.NullLocator())
    
    ax[1].plot(
        np.fft.rfftfreq(len(t), d=np.diff(t)[0]),
        np.abs(np.fft.rfft(theta - np.mean(theta))),
        linewidth=0.5
    )
    ax[1].set_xlabel('$\\nu$ [Hz]')
    ax[1].set_ylabel('$\mathcal F[\\theta](\\nu)$')
    ax[1].set_xlim(0, 1)
    ax[1].yaxis.set_major_locator(plt.NullLocator())
    ax[1].yaxis.set_minor_locator(plt.NullLocator())
    
    fig.savefig(f'generated/pendulum {os.path.basename(filepath)} phase and FT.pdf')
    plt.show(fig)
    plt.close(fig)

data/pendulum/underdamped 300mA 4.4V.csv

data/pendulum/underdamped 300mA 4.5V.csv

data/pendulum/underdamped 300mA 5.5V.csv


def get_H(N, k, a, g):
    H = np.zeros((2*N, 2*N), dtype=complex)
    
    for n in range(N):
        # even layers are connected horizontally in the same unit cell (R=0)
        if n % 2 == 0: H[2*n, 2*n+1] = 1
        # odd layers are connected horizontally periodically in different unit cells (R=a)
        else: H[2*n, 2*n+1] = np.exp(-1j*k*a)
        
        # make the connections vertically between adjacent layers (R=0)
        if n > 0:
            H[2*n-2, 2*n] = 1
            H[2*n-1, 2*n+1] = 1
            
    H = (H + H.conj().T)/2  # make hermitian
    H *= g  # multiply by the hopping integral
    return H


!pip3 install jupyter_contrib_nbextensions ipywidgets
!jupyter nbextension enable --py widgetsnbextension
# Then restart the notebook by stopping
# the process in the terminal (CTRL+C)
# and rerunning `jupyter notebook`

Defaulting to user installation because normal site-packages is not writeable
Requirement already satisfied: jupyter_contrib_nbextensions in /home/nicole/.local/lib/python3.9/site-packages (0.5.1)
Requirement already satisfied: ipywidgets in /home/nicole/.local/lib/python3.9/site-packages (7.6.3)
Requirement already satisfied: widgetsnbextension~=3.5.0 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (3.5.1)
Requirement already satisfied: ipykernel>=4.5.1 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (5.5.0)
Requirement already satisfied: traitlets>=4.3.1 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (5.0.5)
Requirement already satisfied: ipython>=4.0.0 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (7.20.0)
Requirement already satisfied: nbformat>=4.2.0 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (5.1.2)
Requirement already satisfied: jupyterlab-widgets>=1.0.0 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (1.0.0)
Requirement already satisfied: jupyter-nbextensions-configurator>=0.4.0 in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (0.4.1)
Requirement already satisfied: tornado in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (6.1)
Requirement already satisfied: jupyter-contrib-core>=0.3.3 in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (0.3.3)
Requirement already satisfied: nbconvert>=4.2 in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (6.0.7)
Requirement already satisfied: notebook>=4.0 in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (6.2.0)
Requirement already satisfied: traitlets>=4.3.1 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (5.0.5)
Requirement already satisfied: jupyter-core in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (4.7.1)
Requirement already satisfied: ipython-genutils in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (0.2.0)
Requirement already satisfied: pyyaml in /usr/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (5.4.1)
Requirement already satisfied: jupyter-latex-envs>=1.3.8 in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (1.4.6)
Requirement already satisfied: jupyter-highlight-selected-word>=0.1.1 in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (0.2.0)
Requirement already satisfied: lxml in /usr/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (4.6.2)
Requirement already satisfied: tornado in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (6.1)
Requirement already satisfied: jupyter-client in /home/nicole/.local/lib/python3.9/site-packages (from ipykernel>=4.5.1->ipywidgets) (6.1.11)
Requirement already satisfied: traitlets>=4.3.1 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (5.0.5)
Requirement already satisfied: ipython>=4.0.0 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (7.20.0)
Requirement already satisfied: setuptools>=18.5 in /usr/lib/python3.9/site-packages (from ipython>=4.0.0->ipywidgets) (54.1.1)
Requirement already satisfied: jedi>=0.16 in /home/nicole/.local/lib/python3.9/site-packages (from ipython>=4.0.0->ipywidgets) (0.18.0)
Requirement already satisfied: pickleshare in /home/nicole/.local/lib/python3.9/site-packages (from ipython>=4.0.0->ipywidgets) (0.7.5)
Requirement already satisfied: pygments in /usr/lib/python3.9/site-packages (from ipython>=4.0.0->ipywidgets) (2.8.1)
Requirement already satisfied: backcall in /home/nicole/.local/lib/python3.9/site-packages (from ipython>=4.0.0->ipywidgets) (0.2.0)
Requirement already satisfied: traitlets>=4.3.1 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (5.0.5)
Requirement already satisfied: pexpect>4.3 in /home/nicole/.local/lib/python3.9/site-packages (from ipython>=4.0.0->ipywidgets) (4.8.0)
Requirement already satisfied: prompt-toolkit!=3.0.0,!=3.0.1,<3.1.0,>=2.0.0 in /home/nicole/.local/lib/python3.9/site-packages (from ipython>=4.0.0->ipywidgets) (3.0.16)
Requirement already satisfied: decorator in /home/nicole/.local/lib/python3.9/site-packages (from ipython>=4.0.0->ipywidgets) (4.4.2)
Requirement already satisfied: parso<0.9.0,>=0.8.0 in /home/nicole/.local/lib/python3.9/site-packages (from jedi>=0.16->ipython>=4.0.0->ipywidgets) (0.8.1)
Requirement already satisfied: jupyter-core in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (4.7.1)
Requirement already satisfied: python-dateutil>=2.1 in /usr/lib/python3.9/site-packages (from jupyter-client->ipykernel>=4.5.1->ipywidgets) (2.8.1)
Requirement already satisfied: tornado in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (6.1)
Requirement already satisfied: pyzmq>=13 in /home/nicole/.local/lib/python3.9/site-packages (from jupyter-client->ipykernel>=4.5.1->ipywidgets) (22.0.3)
Requirement already satisfied: traitlets>=4.3.1 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (5.0.5)
Requirement already satisfied: notebook>=4.0 in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (6.2.0)
Requirement already satisfied: jupyter-core in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (4.7.1)
Requirement already satisfied: traitlets>=4.3.1 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (5.0.5)
Requirement already satisfied: setuptools>=18.5 in /usr/lib/python3.9/site-packages (from ipython>=4.0.0->ipywidgets) (54.1.1)
Requirement already satisfied: tornado in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (6.1)
Requirement already satisfied: traitlets>=4.3.1 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (5.0.5)
Requirement already satisfied: traitlets>=4.3.1 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (5.0.5)
Requirement already satisfied: notebook>=4.0 in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (6.2.0)
Requirement already satisfied: jupyter-core in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (4.7.1)
Requirement already satisfied: ipython>=4.0.0 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (7.20.0)
Requirement already satisfied: nbconvert>=4.2 in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (6.0.7)
Requirement already satisfied: jupyter-contrib-core>=0.3.3 in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (0.3.3)
Requirement already satisfied: notebook>=4.0 in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (6.2.0)
Requirement already satisfied: traitlets>=4.3.1 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (5.0.5)
Requirement already satisfied: jupyter-core in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (4.7.1)
Requirement already satisfied: pyyaml in /usr/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (5.4.1)
Requirement already satisfied: tornado in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (6.1)
Requirement already satisfied: defusedxml in /home/nicole/.local/lib/python3.9/site-packages (from nbconvert>=4.2->jupyter_contrib_nbextensions) (0.6.0)
Requirement already satisfied: pandocfilters>=1.4.1 in /home/nicole/.local/lib/python3.9/site-packages (from nbconvert>=4.2->jupyter_contrib_nbextensions) (1.4.3)
Requirement already satisfied: testpath in /home/nicole/.local/lib/python3.9/site-packages (from nbconvert>=4.2->jupyter_contrib_nbextensions) (0.4.4)
Requirement already satisfied: jinja2>=2.4 in /home/nicole/.local/lib/python3.9/site-packages (from nbconvert>=4.2->jupyter_contrib_nbextensions) (2.11.3)
Requirement already satisfied: jupyterlab-pygments in /home/nicole/.local/lib/python3.9/site-packages (from nbconvert>=4.2->jupyter_contrib_nbextensions) (0.1.2)
Requirement already satisfied: bleach in /home/nicole/.local/lib/python3.9/site-packages (from nbconvert>=4.2->jupyter_contrib_nbextensions) (3.3.0)
Requirement already satisfied: pygments in /usr/lib/python3.9/site-packages (from ipython>=4.0.0->ipywidgets) (2.8.1)
Requirement already satisfied: jupyter-core in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (4.7.1)
Requirement already satisfied: entrypoints>=0.2.2 in /home/nicole/.local/lib/python3.9/site-packages (from nbconvert>=4.2->jupyter_contrib_nbextensions) (0.3)
Requirement already satisfied: mistune<2,>=0.8.1 in /home/nicole/.local/lib/python3.9/site-packages (from nbconvert>=4.2->jupyter_contrib_nbextensions) (0.8.4)
Requirement already satisfied: nbformat>=4.2.0 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (5.1.2)
Requirement already satisfied: nbclient<0.6.0,>=0.5.0 in /home/nicole/.local/lib/python3.9/site-packages (from nbconvert>=4.2->jupyter_contrib_nbextensions) (0.5.2)
Requirement already satisfied: traitlets>=4.3.1 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (5.0.5)
Requirement already satisfied: six>=1.9.0 in /usr/lib/python3.9/site-packages (from bleach->nbconvert>=4.2->jupyter_contrib_nbextensions) (1.15.0)
Requirement already satisfied: webencodings in /usr/lib/python3.9/site-packages (from bleach->nbconvert>=4.2->jupyter_contrib_nbextensions) (0.5.1)
Requirement already satisfied: packaging in /usr/lib/python3.9/site-packages (from bleach->nbconvert>=4.2->jupyter_contrib_nbextensions) (20.9)
Requirement already satisfied: MarkupSafe>=0.23 in /usr/lib/python3.9/site-packages (from jinja2>=2.4->nbconvert>=4.2->jupyter_contrib_nbextensions) (1.1.1)
Requirement already satisfied: pygments in /usr/lib/python3.9/site-packages (from ipython>=4.0.0->ipywidgets) (2.8.1)
Requirement already satisfied: nbformat>=4.2.0 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (5.1.2)
Requirement already satisfied: jupyter-client in /home/nicole/.local/lib/python3.9/site-packages (from ipykernel>=4.5.1->ipywidgets) (6.1.11)
Requirement already satisfied: async-generator in /home/nicole/.local/lib/python3.9/site-packages (from nbclient<0.6.0,>=0.5.0->nbconvert>=4.2->jupyter_contrib_nbextensions) (1.10)
Requirement already satisfied: nest-asyncio in /home/nicole/.local/lib/python3.9/site-packages (from nbclient<0.6.0,>=0.5.0->nbconvert>=4.2->jupyter_contrib_nbextensions) (1.5.1)
Requirement already satisfied: traitlets>=4.3.1 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (5.0.5)
Requirement already satisfied: traitlets>=4.3.1 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (5.0.5)
Requirement already satisfied: jsonschema!=2.5.0,>=2.4 in /home/nicole/.local/lib/python3.9/site-packages (from nbformat>=4.2.0->ipywidgets) (3.2.0)
Requirement already satisfied: jupyter-core in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (4.7.1)
Requirement already satisfied: ipython-genutils in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (0.2.0)
Requirement already satisfied: six>=1.9.0 in /usr/lib/python3.9/site-packages (from bleach->nbconvert>=4.2->jupyter_contrib_nbextensions) (1.15.0)
Requirement already satisfied: pyrsistent>=0.14.0 in /home/nicole/.local/lib/python3.9/site-packages (from jsonschema!=2.5.0,>=2.4->nbformat>=4.2.0->ipywidgets) (0.17.3)
Requirement already satisfied: attrs>=17.4.0 in /usr/lib/python3.9/site-packages (from jsonschema!=2.5.0,>=2.4->nbformat>=4.2.0->ipywidgets) (20.3.0)
Requirement already satisfied: setuptools>=18.5 in /usr/lib/python3.9/site-packages (from ipython>=4.0.0->ipywidgets) (54.1.1)
Requirement already satisfied: argon2-cffi in /home/nicole/.local/lib/python3.9/site-packages (from notebook>=4.0->jupyter_contrib_nbextensions) (20.1.0)
Requirement already satisfied: Send2Trash>=1.5.0 in /home/nicole/.local/lib/python3.9/site-packages (from notebook>=4.0->jupyter_contrib_nbextensions) (1.5.0)
Requirement already satisfied: jupyter-core in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (4.7.1)
Requirement already satisfied: pyzmq>=13 in /home/nicole/.local/lib/python3.9/site-packages (from jupyter-client->ipykernel>=4.5.1->ipywidgets) (22.0.3)
Requirement already satisfied: traitlets>=4.3.1 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (5.0.5)
Requirement already satisfied: terminado>=0.8.3 in /home/nicole/.local/lib/python3.9/site-packages (from notebook>=4.0->jupyter_contrib_nbextensions) (0.9.2)
Requirement already satisfied: nbformat>=4.2.0 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (5.1.2)
Requirement already satisfied: ipykernel>=4.5.1 in /home/nicole/.local/lib/python3.9/site-packages (from ipywidgets) (5.5.0)
Requirement already satisfied: ipython-genutils in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (0.2.0)
Requirement already satisfied: tornado in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (6.1)
Requirement already satisfied: prometheus-client in /home/nicole/.local/lib/python3.9/site-packages (from notebook>=4.0->jupyter_contrib_nbextensions) (0.9.0)
Requirement already satisfied: jinja2>=2.4 in /home/nicole/.local/lib/python3.9/site-packages (from nbconvert>=4.2->jupyter_contrib_nbextensions) (2.11.3)
Requirement already satisfied: nbconvert>=4.2 in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (6.0.7)
Requirement already satisfied: jupyter-client in /home/nicole/.local/lib/python3.9/site-packages (from ipykernel>=4.5.1->ipywidgets) (6.1.11)
Requirement already satisfied: cffi>=1.0.0 in /usr/lib/python3.9/site-packages (from argon2-cffi->notebook>=4.0->jupyter_contrib_nbextensions) (1.14.5)
Requirement already satisfied: six>=1.9.0 in /usr/lib/python3.9/site-packages (from bleach->nbconvert>=4.2->jupyter_contrib_nbextensions) (1.15.0)
Requirement already satisfied: pycparser in /usr/lib/python3.9/site-packages (from cffi>=1.0.0->argon2-cffi->notebook>=4.0->jupyter_contrib_nbextensions) (2.20)
Requirement already satisfied: pyparsing>=2.0.2 in /usr/lib/python3.9/site-packages (from packaging->bleach->nbconvert>=4.2->jupyter_contrib_nbextensions) (2.4.7)
Requirement already satisfied: ptyprocess>=0.5 in /home/nicole/.local/lib/python3.9/site-packages (from pexpect>4.3->ipython>=4.0.0->ipywidgets) (0.7.0)
Requirement already satisfied: wcwidth in /home/nicole/.local/lib/python3.9/site-packages (from prompt-toolkit!=3.0.0,!=3.0.1,<3.1.0,>=2.0.0->ipython>=4.0.0->ipywidgets) (0.2.5)
Requirement already satisfied: six>=1.9.0 in /usr/lib/python3.9/site-packages (from bleach->nbconvert>=4.2->jupyter_contrib_nbextensions) (1.15.0)
Requirement already satisfied: tornado in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (6.1)
Requirement already satisfied: ptyprocess>=0.5 in /home/nicole/.local/lib/python3.9/site-packages (from pexpect>4.3->ipython>=4.0.0->ipywidgets) (0.7.0)
Requirement already satisfied: ipython-genutils in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (0.2.0)
Requirement already satisfied: notebook>=4.0 in /home/nicole/.local/lib/python3.9/site-packages (from jupyter_contrib_nbextensions) (6.2.0)
Enabling notebook extension jupyter-js-widgets/extension...
      - Validating: OK


def f(x):
    return x


from ipywidgets import interact

interact(f, x=(-10, 10, 1))

<function __main__.f(x)>


def plot_H(k):
    fig, (ax_abs, ax_phase) = plt.subplots(ncols=2, figsize=(8, 3), dpi=150, facecolor='white')

    H = get_H(4, k, 1, -2.7)

    im_abs = ax_abs.matshow(np.abs(H), vmin=-1.5, vmax=1.5)
    ax_abs.set_title('Absolute value')
    fig.colorbar(im_abs, ax=ax_abs)

    im_phase = ax_phase.matshow(np.angle(H), vmin=-np.pi, vmax=np.pi)
    ax_phase.set_title('Phase')
    fig.colorbar(im_phase, ax=ax_phase)

    plt.show(fig)
    plt.close(fig)
    
interact(plot_H, k=(-np.pi, np.pi, np.pi/3))

<function __main__.plot_H(k)>


a = 1
g = -2.7
M = 100

fig, ax = plt.subplots(ncols=4, sharey=True, dpi=150, facecolor='white', figsize=(8, 3))

for i, N in enumerate(range(2, 5+1)):
    k = np.linspace(-np.pi/a, np.pi/a, M)
    E = np.zeros((M, 2*N))  # M samples, 2*N bands
    for j in range(M):
        E[j, :] = np.real(np.linalg.eigvalsh(get_H(N, k[j], a, g)))
    
    ax[i].plot(k, E)
    ax[i].set_xlabel('$k$')
    ax[i].set_ylabel('$E$ [eV]')
    ax[i].set_title(f'$N={N}$')
    
    # pi/a multiple ticks
    ax[i].xaxis.set_major_formatter(
        matplotlib.ticker.FuncFormatter(
           lambda val, pos: '$\\frac{{{:.0g}\pi}}{{a}}$'.format(val/np.pi/a) if val != 0 else '$0$'
        )
    )
    ax[i].xaxis.set_major_locator(matplotlib.ticker.MultipleLocator(base=np.pi/a))

for ax in fig.axes: ax.label_outer()  # Label only outer axes
fig.savefig('generated/bands.pdf')
plt.show(fig)
plt.close(fig)


R = 2
b = 0.4
a0 = 2e4
eps0 = 8.85418782e-12
V0 = 10
N = 100

def eps_r_(r):
    if 0 <= r and r < b: return 1
    else: return 2*(1 + np.arctan(2*np.pi*(r-b)/R))

def rho_lib_(r):
    if 0 <= r and r < b: return eps0*a0*(r**2/b**2*(2*r/b-3)+1)
    else: return 0


eps_r = np.vectorize(eps_r_, otypes=[float])
rho_lib = np.vectorize(rho_lib_, otypes=[float])


fig, axl = plt.subplots(dpi=150, facecolor='white')
axr = axl.twinx()

r = np.linspace(0, R, 1000)

colorl = 'tab:blue'
colorr = 'tab:orange'

axl.plot(r, eps_r(r), color=colorl)
axl.set_xlabel('$r$ [m]')
axl.set_ylabel('$\\epsilon_\\textrm{\\small r}$', color=colorl)
axl.tick_params(axis='y', labelcolor=colorl)

axr.plot(r, rho_lib(r)/(a0*eps0), color=colorr)
axr.set_ylabel('$\\rho_\\textrm{\\small lib}/(a_0 \\epsilon_0)$', color=colorr)
axr.tick_params(axis='y', labelcolor=colorr)

plt.show(fig)
plt.close(fig)


r = np.linspace(0, R, N)
A = np.zeros((N, N))
xi = np.zeros(N)
h = np.diff(r)

# Construct matrix body and rhs
for i in range(len(A)-1):
    xi[i] += h[i]*r[i]*rho_lib(r[i])/2
    xi[i+1] += h[i]*r[i+1]*rho_lib(r[i+1])/2
    A[i, i] += (1.0/h[i])*eps0*(r[i]*eps_r(r[i]) + r[i+1]*eps_r(r[i+1]))/2
    A[i+1, i+1] += (1.0/h[i])*eps0*(r[i]*eps_r(r[i]) + r[i+1]*eps_r(r[i+1]))/2
    A[i, i+1] += (-1.0/h[i])*eps0*(r[i]*eps_r(r[i]) + r[i+1]*eps_r(r[i+1]))/2
    A[i+1, i] += (-1.0/h[i])*eps0*(r[i]*eps_r(r[i]) + r[i+1]*eps_r(r[i+1]))/2

# Boundary conditions
xi[-1] = V0
A[-1, -2] = 0
A[-1, -1] = 1

# Solve
phi = np.linalg.solve(A, xi)


fig, ax = plt.subplots(dpi=150, facecolor='white')

ax.plot(r, phi)

ax.set_xlabel('$r$ [m]')
ax.set_ylabel('$\\phi$ [V]')
fig.savefig('generated/electrostatics potential.pdf')
plt.show(fig)
plt.close(fig)


# The lazy way (incorrect, but is quick to code and less difficult to think about)
# E = -np.gradient(phi, r)
# D = eps_r(r)*eps0*E
# div_E = np.gradient(r*E, r)/r  # recall : cyclindrical coordinates
# div_D = np.gradient(r*D, r)/r
# rho_pol = eps0*div_E - div_D

# The hard way (finite centered differences)
mids = lambda x: (x[1:] + x[:-1]) / 2  # midpoints of an array
diff = lambda x, y: (mids(x), np.diff(y)/np.diff(x))  # derivative
divg = lambda x, y: (mids(x), np.diff(x*y)/np.diff(x)/mids(x))  # cyclindrical divergence
rmid, E = diff(r, -phi)
D = eps_r(rmid)*eps0*E
rmidmid, div_E = divg(rmid, E)
rmidmid, div_D = divg(rmid, D)
rho_pol = eps0*div_E - div_D


fig, ax = plt.subplots(dpi=150, facecolor='white')

ax.plot(rmid, E)

ax.set_xlabel('$r$ [m]')
ax.set_ylabel('$E$ [V/m]')
fig.savefig('generated/electrostatics field.pdf')
plt.show(fig)
plt.close(fig)


fig, ax = plt.subplots(dpi=150, facecolor='white')

ax.plot(rmidmid, rho_lib(rmidmid), label='$\\rho_\\textrm{\\small lib}$ [C/m$^3$]')
ax.plot(rmidmid, rho_pol, label='$\\rho_\\textrm{\\small pol}$ [C/m$^3$]')
ax.plot(rmidmid, div_D, '--', label='$\\vec \\nabla \\cdot \\vec D$ [C/m$^3$]')

ax.set_xlabel('$r$ [m]')
ax.legend(fontsize='small')
fig.savefig('generated/electrostatics charges.pdf')
plt.show(fig)
plt.close(fig)

Jupyter notebooks introduction¶

Commands¶

Markdown¶

Exporting¶

(Optional) theming¶

(Optional) nbextensions¶

(Optional) custom keyboard shortcuts¶

Running this notebook¶

Packages¶

LaTeX¶

Using Google Colab notebooks¶

Magneto-optic effects - Lummer Gehreke interferometer¶

First plots¶

Linear fit¶

Compute Bohr magneton¶

A note on the quality of a fit¶

Residuals and coefficient of determination R2R2R^2¶

Pearson correlation coefficient rrr¶

Display the fit¶

Add the errors !¶

Using uncertainties¶

A plot with errorbars¶

Multiple lines on the same plot¶

Loading all files in a directory which match a pattern¶

Automating the legend cases with lookup dict¶

Using LATEXLATEX\LaTeX¶

Cold plasma - Langmuir probe characteristic curve¶

Log scale¶

Haha, gaussian smoothing go BRRRRR¶

Yeet the hysterisis¶

(Advanced) Regex expressions¶

Chaotic motorized pendulum¶

Boolean masking¶

Using numerical derivatives to sync the signal¶

Fourier transform¶

My first Fourier transform¶

Removing the mean¶

Solid state physics - energy bands of one-dimensional graphene monoribbons¶

Compute the Hamiltonian¶

Interactive visualizations¶

Plot eigenvalues (energies) for different eigenvectors¶

Electrostatics - dielectric cylinder¶

Define and plot ϵrϵr\epsilon_\textrm{r} and ρlibρlib\rho_\textrm{lib}¶

Compute the matrix and solve¶

Compute polarization charges¶

Residuals and coefficient of determination $R^2$ ¶

Pearson correlation coefficient $r$ ¶

Using `uncertainties`¶

Using $\LaTeX$ ¶

Define and plot $\epsilon_\textrm{r}$ and $\rho_\textrm{lib}$ ¶