# Burgers equation

* Physical space 
\begin{align}
u_t + uu_x = \mu u_{xx}
\end{align}
* Spectral space: $\hat{u} = \mathscr{F}\{u\}$
\begin{align}
\hat{u_t} =  -\mu k_x^{2}\hat{u} 
                                       - \mathscr{F}\{ {\mathscr{F}^{-1}\{ \hat{u} \} \mathscr{F}^{-1} \{ i k_x \hat{u} \} } \}
\end{align}


# Imports


```python
import numpy as np
import matplotlib.pyplot as plt
from rkstiff.grids import construct_x_kx_rfft
from rkstiff.derivatives import dx_rfft
from rkstiff.if34 import IF34
from rkstiff.etd35 import ETD35
%matplotlib inline
```

# Helper graph functions


```python
def plotResults(u0,uf,u0FFT,ufFFT):
    """
    Helper function that plots initial 1D-array along with final propagated result.
    Also plots the initial spectral (fft) 1D-array along with final propagated spectral result.
    
    INPUTS
        u0 - initial 1D-array
        uf - final propagated 1D-array
        u0FFT - initial spectral 1D-array
        uFFT - final propagated spectral 1D-array
        
    """
    
    fig = plt.figure(figsize=(12,8))
    ax1 = fig.add_subplot(2,2,1)
    ax2 = fig.add_subplot(2,2,2,sharey = ax1)
    ax3 = fig.add_subplot(2,2,3)
    ax4 = fig.add_subplot(2,2,4,sharey = ax3)
    ax1.plot(x,u0)
    ax2.plot(x,uf)
    ax3.plot(kx,np.abs(u0FFT)**2)
    ax4.plot(kx,np.abs(ufFFT)**2)
    ax3.set_yscale('log')
    ax4.set_yscale('log')
    fig.tight_layout()
    
    
def waterfall(x,t,u,**kwargs):
    if 'figsize' in kwargs:
        fig = plt.figure(figsize=kwargs['figsize'])
    else:
        fig = plt.figure(figsize=(6,6))
    ax = fig.add_subplot(1,1,1,projection='3d')
    ax.w_xaxis.set_pane_color((0,0,0,0))
    ax.w_yaxis.set_pane_color((0,0,0,0))
    ax.w_zaxis.set_pane_color((0,0,0,0))
    for i,snapshot in enumerate(u):
        ax.plot(x,t[i]*np.ones_like(x),snapshot,color='black')
    plt.xlim([x[0],x[-1]])
    plt.ylim([t[0],t[-1]])
    plt.tight_layout()
    return ax
```

# Construct grids


```python
# uniform grid spacing, field assumed to be real-valued -> construct_x_kx_rfft
N = 8192
a,b = -np.pi, np.pi
x,kx = construct_x_kx_rfft(N,a,b)
```

# Linear operator and nonlinear function


```python
mu = 0.0005
L = -mu*kx**2
def NL(uf):
    u = np.fft.irfft(uf)
    ux = np.fft.irfft(1j*kx*uf)
    return -np.fft.rfft(u*ux)
```

# Initial field to be propagated


```python
u0 = np.exp(-10*np.sin(x/2)**2)
u0FFT = np.fft.rfft(u0)
plt.plot(x,u0)
```

# Initialize IF34 solver


```python
solver = IF34(linop=L,NLfunc=NL,epsilon=1e-8)
```

# Propagate step by step


```python
h = 0.01
uFFT = u0FFT.copy()
t = 0
while t < 1:
    uFFT,h,h_suggest = solver.step(uFFT,h)
    t += h
    # use suggested step 
    h = h_suggest  
    
print('tf = ',t)
uf = np.fft.irfft(uFFT)
plotResults(u0,uf,u0FFT,uFFT)
```

# Propagate from time t0 to tf


```python
# store_data -> propagated field stored in solver.u at times solver.t 
# store_freq -> propagated field values stored on every store_freq step (default is every step)
uFFT = solver.evolve(u0FFT,t0=0,tf=1,store_data=True,store_freq=50)
```

# Graph result


```python
U = []
for u in solver.u:
    U.append(np.fft.irfft(u))
U = np.array(U)
t = np.array(solver.t)
ax = waterfall(x,t,U,figsize=(8,8))
ax.grid(False)
ax.axis(False)
ax.view_init(62,-69)
```

d2d2aacdb2d00b5a8f76f903144a5fc023c17c35
ipynb
Jupyter Notebook
demos/burgers.ipynb
whalenpt/rkstiff
9fbec7ddd123cc644d392933b518d342751b4cd8
MIT
2021-11-05T15:35:21.000Z
2022-01-17T10:20:57.000Z
demos/burgers.ipynb
whalenpt/rkstiff
9fbec7ddd123cc644d392933b518d342751b4cd8
MIT
demos/burgers.ipynb
whalenpt/rkstiff
9fbec7ddd123cc644d392933b518d342751b4cd8
MIT
Qwen/Qwen-72B
1. YES
2. YES