# Straight Forward Expansion 


```python
import sympy as sp
from sympy.simplify.fu import TR7, TR8
```


```python
N = 3

# Define the symbolic parameters 
# epsilon = sp.symbols('epsilon_(1:' + str(N+1) + ')')
epsilon = sp.symbols('epsilon')
omega_0 = sp.symbols('omega_0')
alpha = sp.symbols('alpha_(2:' + str(N+1) + ')')

# Define time variable
t = sp.symbols('t', real=True)

# Define function variable 
x = sp.Function('x')(t)
xdot = sp.Derivative(x, t)  # first time derivative 
xddot = sp.Derivative(xdot, t)  # second time derivative

# EOM
EOM = xddot + omega_0**2 * x + sum([alpha[i-2] * x**i for i in range(2,N+1)])
EOM
```




$\displaystyle \alpha_{2} x^{2}{\left(t \right)} + \alpha_{3} x^{3}{\left(t \right)} + \omega_{0}^{2} x{\left(t \right)} + \frac{d^{2}}{d t^{2}} x{\left(t \right)}$




```python
# Exponential expansion of the variable x 
x1 = sp.Function('x_1', real=True)(t)
x2 = sp.Function('x_2', real=True)(t)
x3 = sp.Function('x_3', real=True)(t)
```


```python
x_i = (x1, x2, x3)
x_e = sum([epsilon**i * x_i[i-1] for i in range(1,N+1)])
x_e
```




$\displaystyle \epsilon^{3} \operatorname{x_{3}}{\left(t \right)} + \epsilon^{2} \operatorname{x_{2}}{\left(t \right)} + \epsilon \operatorname{x_{1}}{\left(t \right)}$




```python
# Substitute this into the EOM 
EOM = EOM.subs(x, x_e)
EOM
```




$\displaystyle \alpha_{2} \left(\epsilon^{3} \operatorname{x_{3}}{\left(t \right)} + \epsilon^{2} \operatorname{x_{2}}{\left(t \right)} + \epsilon \operatorname{x_{1}}{\left(t \right)}\right)^{2} + \alpha_{3} \left(\epsilon^{3} \operatorname{x_{3}}{\left(t \right)} + \epsilon^{2} \operatorname{x_{2}}{\left(t \right)} + \epsilon \operatorname{x_{1}}{\left(t \right)}\right)^{3} + \omega_{0}^{2} \left(\epsilon^{3} \operatorname{x_{3}}{\left(t \right)} + \epsilon^{2} \operatorname{x_{2}}{\left(t \right)} + \epsilon \operatorname{x_{1}}{\left(t \right)}\right) + \frac{\partial^{2}}{\partial t^{2}} \left(\epsilon^{3} \operatorname{x_{3}}{\left(t \right)} + \epsilon^{2} \operatorname{x_{2}}{\left(t \right)} + \epsilon \operatorname{x_{1}}{\left(t \right)}\right)$




```python
EOM = sp.expand(sp.expand(EOM).doit())
```


```python
# Collect the coefficients for the epsilons 
epsilon_Eq = sp.collect(EOM, epsilon, evaluate=False)
epsilon_1_Eq = sp.Eq(epsilon_Eq[epsilon], 0)
epsilon_1_Eq
```




$\displaystyle \omega_{0}^{2} \operatorname{x_{1}}{\left(t \right)} + \frac{d^{2}}{d t^{2}} \operatorname{x_{1}}{\left(t \right)} = 0$




```python
epsilon_2_Eq = sp.Eq(epsilon_Eq[epsilon**2], 0)
epsilon_2_Eq
```




$\displaystyle \alpha_{2} \operatorname{x_{1}}^{2}{\left(t \right)} + \omega_{0}^{2} \operatorname{x_{2}}{\left(t \right)} + \frac{d^{2}}{d t^{2}} \operatorname{x_{2}}{\left(t \right)} = 0$




```python
epsilon_3_Eq = sp.Eq(epsilon_Eq[epsilon**3], 0)
epsilon_3_Eq
```




$\displaystyle 2 \alpha_{2} \operatorname{x_{1}}{\left(t \right)} \operatorname{x_{2}}{\left(t \right)} + \alpha_{3} \operatorname{x_{1}}^{3}{\left(t \right)} + \omega_{0}^{2} \operatorname{x_{3}}{\left(t \right)} + \frac{d^{2}}{d t^{2}} \operatorname{x_{3}}{\left(t \right)} = 0$




```python
sp.dsolve(epsilon_1_Eq, x1)
```




$\displaystyle \operatorname{x_{1}}{\left(t \right)} = C_{1} e^{- i \omega_{0} t} + C_{2} e^{i \omega_{0} t}$



Not as convenient as working in the polar form


```python
a = sp.symbols('a')
beta = sp.symbols('beta')
x1_polar = a * sp.cos(omega_0 * t + beta)
x1_polar
```




$\displaystyle a \cos{\left(\beta + \omega_{0} t \right)}$



Updated $\epsilon$-2 equation 


```python
epsilon_2_Eq = sp.expand(TR7(epsilon_2_Eq.subs(x1, x1_polar)))
epsilon_2_Eq
```




$\displaystyle \frac{a^{2} \alpha_{2} \cos{\left(2 \beta + 2 \omega_{0} t \right)}}{2} + \frac{a^{2} \alpha_{2}}{2} + \omega_{0}^{2} \operatorname{x_{2}}{\left(t \right)} + \frac{d^{2}}{d t^{2}} \operatorname{x_{2}}{\left(t \right)} = 0$




```python
sp.dsolve(epsilon_2_Eq, x2)
```




$\displaystyle \operatorname{x_{2}}{\left(t \right)} = C_{1} e^{- i \omega_{0} t} + C_{2} e^{i \omega_{0} t} + \frac{a^{2} \alpha_{2} \cos{\left(2 \beta + 2 \omega_{0} t \right)}}{6 \omega_{0}^{2}} - \frac{a^{2} \alpha_{2}}{2 \omega_{0}^{2}}$



We only want to keep the particular solution in the above


```python
x2_p = alpha[0] * a**2 * (sp.cos(2 * omega_0 * t + 2*beta) - 3)/6/omega_0**2
x2_p
```




$\displaystyle \frac{a^{2} \alpha_{2} \left(\cos{\left(2 \beta + 2 \omega_{0} t \right)} - 3\right)}{6 \omega_{0}^{2}}$



Updated $\epsilon$-3 equation


```python
epsilon_3_Eq = TR7(epsilon_3_Eq.subs([
    (x1, x1_polar), (x2, x2_p)
]))
epsilon_3_Eq = TR8(epsilon_3_Eq)
epsilon_3_Eq
```




$\displaystyle \frac{a^{3} \alpha_{2}^{2} \left(- 5 \cos{\left(\beta + \omega_{0} t \right)} + \cos{\left(3 \beta + 3 \omega_{0} t \right)}\right)}{6 \omega_{0}^{2}} + a^{3} \alpha_{3} \left(\frac{3 \cos{\left(\beta + \omega_{0} t \right)}}{4} + \frac{\cos{\left(3 \beta + 3 \omega_{0} t \right)}}{4}\right) + \omega_{0}^{2} \operatorname{x_{3}}{\left(t \right)} + \frac{d^{2}}{d t^{2}} \operatorname{x_{3}}{\left(t \right)} = 0$




```python
sp.dsolve(epsilon_3_Eq)
```




$\displaystyle \operatorname{x_{3}}{\left(t \right)} = C_{1} e^{- i \omega_{0} t} + C_{2} e^{i \omega_{0} t} + \frac{5 a^{3} \alpha_{2}^{2} t \sin{\left(\beta + \omega_{0} t \right)}}{12 \omega_{0}^{3}} + \frac{a^{3} \alpha_{2}^{2} \cos{\left(3 \beta + 3 \omega_{0} t \right)}}{48 \omega_{0}^{4}} - \frac{3 a^{3} \alpha_{3} t \sin{\left(\beta + \omega_{0} t \right)}}{8 \omega_{0}} + \frac{a^{3} \alpha_{3} \cos{\left(3 \beta + 3 \omega_{0} t \right)}}{32 \omega_{0}^{2}}$



fed8b6634941b0ef171636c033a14452d8f0500f
ipynb
Jupyter Notebook
notebooks/straight_forward_expansion.ipynb
smallpondtom/pNLsys
ebaea9c2945391fccd076f8006baee7066951097
MIT
2021-12-17T16:44:14.000Z
2022-02-22T05:40:43.000Z
notebooks/straight_forward_expansion.ipynb
smallpondtom/pNLsys
ebaea9c2945391fccd076f8006baee7066951097
MIT
notebooks/straight_forward_expansion.ipynb
smallpondtom/pNLsys
ebaea9c2945391fccd076f8006baee7066951097
MIT
Qwen/Qwen-72B
1. YES
2. YES