# Automatic Transformation of Nth Order ODEs to 1st Order ODEs

ModelingToolkit has a system for transformations of mathematical systems. These transformations allow for symbolically changing the representation of the model to problems that are easier to numerically solve. One simple to demonstrate transformation is the `ode_order_lowering`

transformation that sends an Nth order ODE to a 1st order ODE.

To see this, let's define a second order riff on the Lorenz equations. We utilize the derivative operator twice here to define the second order:

```
using ModelingToolkit, OrdinaryDiffEq
@parameters σ ρ β
@variables t x(t) y(t) z(t)
D = Differential(t)
eqs = [D(D(x)) ~ σ*(y-x),
D(y) ~ x*(ρ-z)-y,
D(z) ~ x*y - β*z]
@named sys = ODESystem(eqs)
```

\[ \begin{align} \mathrm{\frac{d}{d t}}\left( \frac{dx(t)}{dt} \right) =& \sigma \left( - x\left( t \right) + y\left( t \right) \right) \\ \frac{dy(t)}{dt} =& - y\left( t \right) + \left( \rho - z\left( t \right) \right) x\left( t \right) \\ \frac{dz(t)}{dt} =& x\left( t \right) y\left( t \right) - \beta z\left( t \right) \end{align} \]

Note that we could've used an alternative syntax for 2nd order, i.e. `D = Differential(t)^2`

and then `E(x)`

would be the second derivative, and this syntax extends to `N`

-th order. Also, we can use `*`

or `∘`

to compose `Differential`

s, like `Differential(t) * Differential(x)`

.

Now let's transform this into the `ODESystem`

of first order components. We do this by simply calling `ode_order_lowering`

:

`sys = ode_order_lowering(sys)`

\[ \begin{align} \frac{dxˍt(t)}{dt} =& \sigma \left( - x\left( t \right) + y\left( t \right) \right) \\ \frac{dy(t)}{dt} =& - y\left( t \right) + \left( \rho - z\left( t \right) \right) x\left( t \right) \\ \frac{dz(t)}{dt} =& x\left( t \right) y\left( t \right) - \beta z\left( t \right) \\ \frac{dx(t)}{dt} =& xˍt\left( t \right) \end{align} \]

Now we can directly numerically solve the lowered system. Note that, following the original problem, the solution requires knowing the initial condition for `x'`

, and thus we include that in our input specification:

```
u0 = [D(x) => 2.0,
x => 1.0,
y => 0.0,
z => 0.0]
p = [σ => 28.0,
ρ => 10.0,
β => 8/3]
tspan = (0.0,100.0)
prob = ODEProblem(sys,u0,tspan,p,jac=true)
sol = solve(prob,Tsit5())
using Plots; plot(sol,vars=(x,y))
```