# RLC Circuit Using Components

So far, we've created models strictly from equations. Now we are going to change gears a bit and create the same RLC circuit but this time using acausal component models.

As with the other examples, there a few bits we always need to do, i.e.,

using ModelingToolkit
using DifferentialEquations
using Plots

In addition, we need to add one more package called Unitful (which you'll need to add as a dependency) for assigning units to variables and values...

using Unitful

And, of course, we define the time, t, and the differentiation operator, D.

@variables t [unit = u"s"]
D = Differential(t)

Note that these will eventually be moving into a package so you won't need to define them for every model. But, for now, we still need to define them.

As in Modelica, ModelingToolkit uses connectors to define the interactions between components. These need to be defined once for each domain and then can be reused. Connectors are the foundation for modeling a given domain. For an electrical system, the connectors are defined as follows:

@connector Pin begin
    v(t), [unit = u"V"]
    i(t), [connect = Flow, unit = u"A"]
end

The @connector macro is used to define a connector (in this case named Pin). Within the connector we have two variables v(t) and i(t). Both of these have units associated with them and the current, i(t), is also marked with connect = Flow. This is analogous to the flow keyword in Modelica.

With the connector definition completed, we've laid the foundation for creating electrical components.

With the connectors defined, we can go through and implement all the components we require.

Our Ground model will be the simplest model we create because it requires only a single equation:

@mtkmodel Ground begin
    @components begin
        g = Pin()
    end
    @equations begin
        g.v ~ 0
    end
end

Because this is the definition of a component model, it uses the @mtkmodel macros (unlike the connectors which start with @connector). The @components section lists all the subcomponents of this model. In this case, there is only a single Pin instance named g. In the equation section, we have a single equation which simply pins the voltage to zero.

Our Resistor model is slightly more sophisticated than our Ground model and it adds a few more complications:

@mtkmodel Resistor begin
    @components begin
        p = Pin()
        n = Pin()
    end
    @variables begin
        v(t), [unit = u"V"]
        i(t), [unit = u"A"]
    end
    @parameters begin
        R, [unit = u"Ω"]
    end
    @equations begin
        v ~ p.v - n.v
        0 ~ p.i + n.i
        i ~ p.i
        v ~ i * R
    end
end

One complication is that it has two subcomponents, the two pins p and n. In addition, it has local variables v(t) and i(t). Finally, in addition to having multiple equations, this component also has parameters defined in the @parameters section.

Our Capacitor model is quite similar to our Resistor model. This will become important in the next section. But for now, just notice that the Capacitor model is pretty much as you would expect given the way our Resistor is defined:

@mtkmodel Capacitor begin
    @components begin
        p = Pin()
        n = Pin()
    end
    @variables begin
        v(t), [unit = u"V"]
        i(t), [unit = u"A"]
    end
    @parameters begin
        C, [unit = u"F"]
    end
    @equations begin
        v ~ p.v - n.v
        0 ~ p.i + n.i
        i ~ p.i
        C * D(v) ~ i
    end
end

Now from a software engineering perspective, there is something very bad about our previous Resistor and Capacitor models. Specifically, there is way too much duplicate code. One of the key principles in software development is the DRY principle which stands for "Don't repeat yourself".

Fortunately, MTK has a way for us to avoid repeating ourselves. We do this by defining a special model which we will call TwoPin that represents all the common bits between Resistor and Capacitor (and, as we will see shortly, Inductor).

@mtkmodel TwoPin begin
    @components begin
        p = Pin()
        n = Pin()
    end
    @variables begin
        v(t), [unit = u"V"]
        i(t), [unit = u"A"]
    end
    @equations begin
        v ~ p.v - n.v
        0 ~ p.i + n.i
        i ~ p.i
    end
end

Note how the @components, @variables and @equations sections here contain the common bits from Resistor and Capacitor. In subsequent sections we will see how this makes the definition of our component models much simpler.

Now that we've defined TwoPin, lets make another pass at our Resistor model:

@mtkmodel Resistor begin
    @extend v, i = twopin = TwoPin()
    @parameters begin
        R, [unit = u"Ω"]
    end
    @equations begin
        v ~ i * R
    end
end

Note how much simpler the model is in this case. We use the @extend macro to extract the variables v and i from the TwoPin instance so we can reference them locally in our equations. All the other sections of TwoPin get merged into the sections in Resistor so that or Resistor also has the connectors and equations of TwoPin.

We can give our Capacitor model the same treatment:

@mtkmodel Capacitor begin
    @extend v, i = twopin = TwoPin()
    @parameters begin
        C, [unit = u"F"]
    end
    @equations begin
        C * D(v) ~ i
    end
end

From a software engineering point, the essential point here is that this model is primarily composed of things that are unique to the capacitor model and avoids repeating anything that it shares with other models.

By leveraging the TwoPin model, we are able to make an Inductor model just as succinctly as we did our latest Resistor and Capacitor models:

@mtkmodel Inductor begin
    @extend v, i = twopin = TwoPin()
    @parameters begin
        L, [unit = u"H"]
    end
    @equations begin
        L * D(i) ~ v
    end
end

In order to complete a component oriented version of the model we created from equations, we need a VoltageSource:

@mtkmodel ConstantVoltage begin
    @extend (v,) = twopin = TwoPin()
    @parameters begin
        V, [unit = u"V"]
    end
    @equations begin
        V ~ v
    end
end

Just as with our previous models, we leverage the TwoPin model. But, in this case, we don't actually need the current variable from TwoPin because the behavior of the current source doesn't depend on current (unlike all the other components previous derived from TwoPin)

With our component models all defined, we can finally assemble our system level model. It uses the same set of macros (e.g., @mtkmodel, @equations, @components) as our other models:

@mtkmodel RLCModel begin
    @components begin
        resistor = Resistor(R = 100.0)
        capacitor = Capacitor(C = 1.0e-3)
        inductor = Inductor(L = 1.0)
        source = ConstantVoltage(V = 24.0)
        ground = Ground()
    end
    @equations begin
        connect(source.p, inductor.p)
        connect(inductor.n, resistor.p)
        connect(inductor.n, capacitor.p)
        connect(resistor.n, ground.g)
        connect(capacitor.n, ground.g)
        connect(source.n, ground.g)
    end
end

With our system model defined, a few more lines are all that is required to simulate the system:

@mtkbuild rlc_model = RLCModel()
u0 = [
    rlc_model.capacitor.v => 0.0,
    rlc_model.inductor.i => 0.0
]
prob = ODEProblem(rlc_model, u0, (0, 1.0))
sol = solve(prob)
display(plot(sol, idxs=[rlc_model.capacitor.v, rlc_model.inductor.i]))

The results from our component oriented model:

...agree exactly with our equation oriented model:

Granted, we had to write a lot more code for the component oriented models. But the difference here is that we didn't have to apply all the conservation laws ourselves and then formulate the resulting ordinary differential equations. All that work was taken care of here. This means that we can combine the models above and any others that we chose to create and not have to worry about errors resulting from improperly formulating the conservation equations. These conservation equations can be derived automatically from the connect equations in our model.

What this means in practice is that our system results in the following equations:

\begin{align}
capacitor_{+}p_{+}v\left( t \right) =& capacitor_{+}v\left( t \right) \\
source_{+}v\left( t \right) =& source_{+}V \\
inductor_{+}p_{+}i\left( t \right) =& inductor_{+}i\left( t \right) \\
ground_{+}g_{+}v\left( t \right) =& 0 \\
capacitor_{+}n_{+}v\left( t \right) =& 0 \\
source_{+}n_{+}v\left( t \right) =& 0 \\
resistor_{+}n_{+}v\left( t \right) =& 0 \\
resistor_{+}v\left( t \right) =& capacitor_{+}p_{+}v\left( t \right) \\
inductor_{+}v\left( t \right) =&  - capacitor_{+}v\left( t \right) + source_{+}v\left( t \right) \\
source_{+}p_{+}v\left( t \right) =& source_{+}v\left( t \right) \\
inductor_{+}p_{+}v\left( t \right) =& source_{+}v\left( t \right) \\
inductor_{+}n_{+}i\left( t \right) =&  - inductor_{+}p_{+}i\left( t \right) \\
source_{+}i\left( t \right) =&  - inductor_{+}p_{+}i\left( t \right) \\
resistor_{+}p_{+}v\left( t \right) =& resistor_{+}v\left( t \right) \\
resistor_{+}i\left( t \right) =& \frac{resistor_{+}v\left( t \right)}{resistor_{+}R} \\
inductor_{+}n_{+}v\left( t \right) =& resistor_{+}v\left( t \right) \\
source_{+}n_{+}i\left( t \right) =&  - source_{+}i\left( t \right) \\
resistor_{+}n_{+}i\left( t \right) =&  - resistor_{+}i\left( t \right) \\
capacitor_{+}i\left( t \right) =& inductor_{+}i\left( t \right) - resistor_{+}i\left( t \right) \\
source_{+}p_{+}i\left( t \right) =&  - source_{+}n_{+}i\left( t \right) \\
resistor_{+}p_{+}i\left( t \right) =&  - resistor_{+}n_{+}i\left( t \right) \\
capacitor_{+}n_{+}i\left( t \right) =&  - capacitor_{+}i\left( t \right) \\
ground_{+}g_{+}i\left( t \right) =&  - inductor_{+}p_{+}i\left( t \right) + resistor_{+}i\left( t \right) + capacitor_{+}i\left( t \right) \\
capacitor_{+}p_{+}i\left( t \right) =&  - capacitor_{+}n_{+}i\left( t \right)
\end{align}

...but MTK automatically reduces this down to:

\begin{align}
\frac{\mathrm{d} capacitor_{+}v\left( t \right)}{\mathrm{d}t} =& \frac{capacitor_{+}i\left( t \right)}{capacitor_{+}C} \\
\frac{\mathrm{d} inductor_{+}i\left( t \right)}{\mathrm{d}t} =& \frac{inductor_{+}v\left( t \right)}{inductor_{+}L}
\end{align}

So we get the best of both worlds. The process of building the models removes many of the tedious, time consuming and error prone aspects of model development thanks to our component oriented approach but MTK's symbolic processing means we get an optimized system of equations to solve.

Note, the complete system of equations for our system can be extracted by executing latexify(full_equations(rlc_model)) and the reduce system of equations can be extracted with latexify(rlc_model), in case you were curious.

Note, much of what we have defined here is already defined in the MTK standard library. But the goal of this example was to demonstrate creation of a simple electrical component library from scratch to drive home the point that these component models are not somehow "built in" to MTK but that anybody can create such components for any engineering domain.

using ModelingToolkit
using DifferentialEquations
using Plots
using Unitful

@variables t [unit = u"s"]
D = Differential(t)

@connector Pin begin
    v(t), [unit = u"V"]
    i(t), [connect = Flow, unit = u"A"]
end

@mtkmodel Ground begin
    @components begin
        g = Pin()
    end
    @equations begin
        g.v ~ 0
    end
end

@mtkmodel TwoPin begin
    @components begin
        p = Pin()
        n = Pin()
    end
    @variables begin
        v(t), [unit = u"V"]
        i(t), [unit = u"A"]
    end
    @equations begin
        v ~ p.v - n.v
        0 ~ p.i + n.i
        i ~ p.i
    end
end

@mtkmodel Resistor begin
    @extend v, i = twopin = TwoPin()
    @parameters begin
        R, [unit = u"Ω"]
    end
    @equations begin
        v ~ i * R
    end
end

@mtkmodel Capacitor begin
    @extend v, i = twopin = TwoPin()
    @parameters begin
        C, [unit = u"F"]
    end
    @equations begin
        C * D(v) ~ i
    end
end

@mtkmodel Inductor begin
    @extend v, i = twopin = TwoPin()
    @parameters begin
        L, [unit = u"H"]
    end
    @equations begin
        L * D(i) ~ v
    end
end

@mtkmodel ConstantVoltage begin
    @extend (v,) = twopin = TwoPin()
    @parameters begin
        V, [unit = u"V"]
    end
    @equations begin
        V ~ v
    end
end

@mtkmodel RLCModel begin
    @components begin
        resistor = Resistor(R = 100.0)
        capacitor = Capacitor(C = 1.0e-3)
        inductor = Inductor(L = 1.0)
        source = ConstantVoltage(V = 24.0)
        ground = Ground()
    end
    @equations begin
        connect(source.p, inductor.p)
        connect(inductor.n, resistor.p)
        connect(inductor.n, capacitor.p)
        connect(resistor.n, ground.g)
        connect(capacitor.n, ground.g)
        connect(source.n, ground.g)
    end
end

@mtkbuild rlc_model = RLCModel()
u0 = [
    rlc_model.capacitor.v => 0.0,
    rlc_model.inductor.i => 0.0
]
prob = ODEProblem(rlc_model, u0, (0, 1.0))
sol = solve(prob)
display(plot(sol, idxs=[rlc_model.capacitor.v, rlc_model.inductor.i]))