# coding: utf-8 # In[1]: # Custom styling, borrowed from https://github.com/ellisonbg/talk-2013-scipy from IPython import display # Make plots centered. H/T to http://stackoverflow.com/a/27168595/1068170 STYLE = """ """ SLIDES_STYLE = """ """ STYLE = STYLE + '\n' + SLIDES_STYLE display.display(display.HTML(STYLE)) # $$u_t + \left[f(u)\right]_x = 0$$ #

# Finite Volume

# $$w(x, t) = \frac{1}{\Delta x} \int_{x - \Delta x / 2}^{x + \Delta x / 2} u(s, t) \, ds$$ # \begin{align*} # w_t &= \frac{1}{\Delta x} \int_{x - \Delta x / 2}^{x + \Delta x / 2} u_t \, ds \\ # &= \frac{1}{\Delta x} \int_{x - \Delta x / 2}^{x + \Delta x / 2} -\left[f(u)\right]_s \, ds \\ # &= -\frac{f\left(u\left(x + \Delta x / 2, t\right)\right) - f\left(u\left(x - \Delta x / 2, t\right)\right)}{\Delta x} # \end{align*} # \begin{align*} # u_j^n &= u\left(j \Delta x, n \Delta t\right) \\ # w_j^n &= w\left(j \Delta x, n \Delta t\right) # \end{align*} # $$\frac{d}{dt} w_j^n = -\frac{1}{\Delta x}\left(f\left(u_{j + \frac{1}{2}}\right) - f\left(u_{j - \frac{1}{2}}\right)\right)$$ # $$\left\{w_j^n\right\}_j \overset{?}{\longrightarrow} \left\{u_{j + \frac{1}{2}}^n\right\}_j$$ # > WENO interpolation ... [is] used ... to transfer information from one domain to another in a high order, nonoscillatory fashion # # ![WENO Transfer Quote](weno_transfer_quote.png "WENO Paper Quote") #

# Finite Difference

# $$\frac{d}{dt} u_j^n = -\frac{1}{\Delta x}\left(\widehat{f}_{j + \frac{1}{2}} - \widehat{f}_{j - \frac{1}{2}}\right)$$ # Define $h(x)$ implicitly via # $$\frac{1}{\Delta x} \int_{x - \Delta x / 2}^{x + \Delta x / 2} h(s) \, ds = f(u(x))$$ # $$\frac{1}{\Delta x} \int_{x - \Delta x / 2}^{x + \Delta x / 2} h(s) \, ds = f(u(x))$$ # $$\Longrightarrow \frac{1}{\Delta x} \left(h\left(x_{j + \frac{1}{2}}\right) -h\left(x_{j - \frac{1}{2}}\right)\right) = \left[f(u)\right]_x$$ # $$\Longrightarrow \frac{1}{\Delta x} \left(h\left(x_{j + \frac{1}{2}}\right) -h\left(x_{j - \frac{1}{2}}\right)\right) \approx \frac{1}{\Delta x}\left(\widehat{f}_{j + \frac{1}{2}} - \widehat{f}_{j - \frac{1}{2}}\right)$$ # $$\widehat{f}_{j + \frac{1}{2}} = h\left(x_{j + \frac{1}{2}}\right) = h_{j + \frac{1}{2}}$$ # $$\overline{h}_j^n = \frac{1}{\Delta x} \int_{x - \Delta x / 2}^{x + \Delta x / 2} h(s) \, ds = f\left(u_j^n\right)$$ # $$\left\{\overline{h}_j^n\right\}_j \overset{?}{\longrightarrow} \left\{h_{j + \frac{1}{2}}^n\right\}_j$$ # In[2]: get_ipython().run_line_magic('matplotlib', 'inline') import seaborn import sympy import weno_computations seaborn.set_palette('husl') # In[3]: weno_computations.make_intro_plots(0) # In[4]: weno_computations.make_intro_plots(1) # In[5]: weno_computations.make_intro_plots(2) # In[6]: weno_computations.make_intro_plots(3) # In[7]: weno_computations.make_intro_plots(4) # In[8]: weno_computations.make_intro_plots(5) # In[9]: weno_computations.make_intro_plots(6) # In[10]: weno_computations.make_intro_plots(7) # In[11]: weno_computations.make_intro_plots(8) # In[12]: weno_computations.make_intro_plots(9) # In[13]: weno_computations.make_intro_plots(10) # In[14]: weno_computations.make_intro_plots(11) # In[15]: weno_computations.make_intro_plots(12) # In[16]: # The slide below can be produced with the following: (approx_minus2, approx_minus1, approx_zero, approx_all) = weno_computations.interp_simple_stencils() approx_minus2 = approx_minus2.replace('=', '&=') approx_minus1 = approx_minus1.replace('=', '&=') approx_zero = approx_zero.replace('=', '&=') approx_all = approx_all.replace('=', '&=') all_math = r"""\begin{align*} %s \\ %s \\ %s \\ %s \end{align*} """ % (approx_minus2, approx_minus1, approx_zero, approx_all) # $$u^{{(1)}}_{{j + \frac{1}{2}}} = \frac{1}{3}\overline{u}_{j-2} - \frac{7}{6}\overline{u}_{j-1} + \frac{11}{6}\overline{u}_{j}$$ # $$u^{{(2)}}_{{j + \frac{1}{2}}} = -\frac{1}{6}\overline{u}_{j-1} + \frac{5}{6}\overline{u}_{j} + \frac{1}{3}\overline{u}_{j+1}$$ # $$u^{{(3)}}_{{j + \frac{1}{2}}} = \frac{1}{3}\overline{u}_{j} + \frac{5}{6}\overline{u}_{j+1} - \frac{1}{6}\overline{u}_{j+2}$$ # $$u^{{(\ast)}}_{{j + \frac{1}{2}}} = \frac{1}{30}\overline{u}_{j-2} - \frac{13}{60}\overline{u}_{j-1} + \frac{47}{60}\overline{u}_{j} + \frac{9}{20}\overline{u}_{j+1} -\frac{1}{20}\overline{u}_{j+2}$$ # $$\frac{1}{3}\overline{u}_{j-2} - \frac{7}{6}\overline{u}_{j-1} + \frac{11}{6}\overline{u}_{j} = u\left(x_{j + \frac{1}{2}}\right) - \frac{\Delta x^3 u_0'''}{4} + \mathcal{O}\left(\Delta x^4\right)$$ # $$-\frac{1}{6}\overline{u}_{j-1} + \frac{5}{6}\overline{u}_{j} + \frac{1}{3}\overline{u}_{j+1} = u\left(x_{j + \frac{1}{2}}\right) + \frac{\Delta x^3 u_0'''}{12} + \mathcal{O}\left(\Delta x^4\right)$$ # $$\frac{1}{3}\overline{u}_{j} + \frac{5}{6}\overline{u}_{j+1} - \frac{1}{6}\overline{u}_{j+2} = u\left(x_{j + \frac{1}{2}}\right) - \frac{\Delta x^3 u_0'''}{12} + \mathcal{O}\left(\Delta x^4\right)$$ # $$u^{{(\ast)}}_{{j + \frac{1}{2}}} = u\left(x_{j + \frac{1}{2}}\right) - \frac{\Delta x^5 u_0'''}{60} + \mathcal{O}\left(\Delta x^6\right)$$ # $$\left(\Delta x\right) \overline{u} = \int_a^b u(x) \, dx = u_0(b - a) + u_0' \left(\frac{b^2 - a^2}{2}\right) + u_0'' \left(\frac{b^3 - a^3}{6}\right) + \cdots$$ # Putting the **W** in **WENO**: # $$u^{{(\ast)}}_{{j + \frac{1}{2}}} = \frac{1}{10} u^{{(1)}}_{{j + \frac{1}{2}}} + \frac{6}{10} u^{{(2)}}_{{j + \frac{1}{2}}} + \frac{3}{10} u^{{(3)}}_{{j + \frac{1}{2}}}$$ # What about the other letters? # What about shocks? # In[17]: weno_computations.discontinuity_to_volume_single_cell(0) # In[18]: weno_computations.discontinuity_to_volume_single_cell(1) # In[19]: weno_computations.make_shock_plot_single_cell() # $$\frac{1}{10} u^{{(1)}}_{{j + \frac{1}{2}}} + \frac{6}{10} u^{{(2)}}_{{j + \frac{1}{2}}} + \frac{3}{10} u^{{(3)}}_{{j + \frac{1}{2}}} \longrightarrow \frac{2}{3} u^{{(2)}}_{{j + \frac{1}{2}}} + \frac{1}{3} u^{{(3)}}_{{j + \frac{1}{2}}}$$ # In[20]: weno_computations.discontinuity_to_volume() # In[21]: weno_computations.make_shock_plot() # Smoothness indicator when approximating a point in the interval $I = \left(x_j - \frac{\Delta x}{2}, x_j + \frac{\Delta x}{2}\right]$ # $$\beta = \sum_k \Delta x^{2k - 1} \int_I \left(\frac{d^k}{dx^k} p\right)^2 dx$$ # Large $\beta$ indicates that the approximation $p(x)$ is not smooth nearby our point. # Recall the **W** in **WENO**: # $$u^{{(\ast)}}_{{j + \frac{1}{2}}} = \frac{1}{10} u^{{(1)}}_{{j + \frac{1}{2}}} + \frac{6}{10} u^{{(2)}}_{{j + \frac{1}{2}}} + \frac{3}{10} u^{{(3)}}_{{j + \frac{1}{2}}}$$ # Rather than using these directly, incorporate a penalty for each polynomial approximation's (lack of) smoothness: $\beta^{(1)}$, $\beta^{(2)}$, $\beta^{(3)}$. # $$\widetilde{w}^{(1)} = \frac{0.1}{\left(\varepsilon + \beta^{(1)}\right)^2}, \; # \widetilde{w}^{(2)} = \frac{0.6}{\left(\varepsilon + \beta^{(2)}\right)^2}, \; # \widetilde{w}^{(3)} = \frac{0.3}{\left(\varepsilon + \beta^{(3)}\right)^2}$$ # $$w^{(j)} = \frac{\widetilde{w}^{(j)}}{\widetilde{w}^{(1)} + \widetilde{w}^{(2)} + \widetilde{w}^{(3)}}$$ # $$w^{(1)} u^{{(1)}}_{{j + \frac{1}{2}}} + w^{(2)} u^{{(2)}}_{{j + \frac{1}{2}}} + w^{(3)} u^{{(3)}}_{{j + \frac{1}{2}}}$$ #

# Back to Conservation Law

# $$\dot{\overline{u}} = -\frac{1}{\Delta x}\left(\widehat{f}\left(u_{j + \frac{1}{2}}^-, u_{j + \frac{1}{2}}^+\right) - \widehat{f}\left(u_{j - \frac{1}{2}}^-, u_{j - \frac{1}{2}}^+\right)\right)$$ # $$u_{j + \frac{1}{2}}^- \longrightarrow \text{Reconstruct from } \overline{u}_{j - 2}, \overline{u}_{j - 1}, \overline{u}_{j}, \overline{u}_{j + 1}, \overline{u}_{j + 2}$$ # $$u_{j + \frac{1}{2}}^+ \longrightarrow \text{Reconstruct from } \overline{u}_{j - 1}, \overline{u}_{j}, \overline{u}_{j + 1}, \overline{u}_{j + 2}, \overline{u}_{j + 3}$$ # $$\dot{\overline{u}} = L\left(\overline{u}\right)$$ # TV-diminishing / strong stability preserving RK scheme # # $$\begin{array}{c | c c c} # 0 & & & \\ # 1 & 1 & & \\ # 1/2 & 1/4 & 1/4 & \\ # \hline # & 1/6 & 1/6 & 2/3 # \end{array}$$ # ![Stability Preserving](tvd_scheme.png "TVD RK Scheme") #

# Example

# ![WENO Example Part 1](traffic_flow_example1.png "WENO Traffic Flow") # ![WENO Example Part 2](traffic_flow_example2.png "WENO Traffic Flow")