-
-
Notifications
You must be signed in to change notification settings - Fork 213
/
Copy pathindex.html
245 lines (227 loc) · 24.8 KB
/
index.html
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
<!DOCTYPE html>
<html lang="en"><head><meta charset="UTF-8"/><meta name="viewport" content="width=device-width, initial-scale=1.0"/><title>Acausal Component-Based Modeling the RC Circuit · ModelingToolkit.jl</title><script data-outdated-warner src="../../assets/warner.js"></script><link rel="canonical" href="https://mtk.sciml.ai/stable/tutorials/acausal_components/"/><link href="https://cdnjs.cloudflare.com/ajax/libs/lato-font/3.0.0/css/lato-font.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/juliamono/0.039/juliamono-regular.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.3/css/fontawesome.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.3/css/solid.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.3/css/brands.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/KaTeX/0.13.11/katex.min.css" rel="stylesheet" type="text/css"/><script>documenterBaseURL="../.."</script><script src="https://cdnjs.cloudflare.com/ajax/libs/require.js/2.3.6/require.min.js" data-main="../../assets/documenter.js"></script><script src="../../siteinfo.js"></script><script src="../../../versions.js"></script><link class="docs-theme-link" rel="stylesheet" type="text/css" href="../../assets/themes/documenter-dark.css" data-theme-name="documenter-dark" data-theme-primary-dark/><link class="docs-theme-link" rel="stylesheet" type="text/css" href="../../assets/themes/documenter-light.css" data-theme-name="documenter-light" data-theme-primary/><script src="../../assets/themeswap.js"></script><link href="../../assets/favicon.ico" rel="icon" type="image/x-icon"/></head><body><div id="documenter"><nav class="docs-sidebar"><a class="docs-logo" href="../../"><img src="../../assets/logo.png" alt="ModelingToolkit.jl logo"/></a><div class="docs-package-name"><span class="docs-autofit"><a href="../../">ModelingToolkit.jl</a></span></div><form class="docs-search" action="../../search/"><input class="docs-search-query" id="documenter-search-query" name="q" type="text" placeholder="Search docs"/></form><ul class="docs-menu"><li><a class="tocitem" href="../../">Home</a></li><li><span class="tocitem">Symbolic Modeling Tutorials</span><ul><li><a class="tocitem" href="../ode_modeling/">Composing Ordinary Differential Equations</a></li><li class="is-active"><a class="tocitem" href>Acausal Component-Based Modeling the RC Circuit</a><ul class="internal"><li><a class="tocitem" href="#Copy-Paste-Example"><span>Copy-Paste Example</span></a></li><li><a class="tocitem" href="#Explanation"><span>Explanation</span></a></li><li><a class="tocitem" href="#Simplifying-and-Solving-this-System"><span>Simplifying and Solving this System</span></a></li></ul></li><li><a class="tocitem" href="../higher_order/">Automatic Transformation of Nth Order ODEs to 1st Order ODEs</a></li><li><a class="tocitem" href="../tearing_parallelism/">Exposing More Parallelism By Tearing Algebraic Equations in ODESystems</a></li><li><a class="tocitem" href="../nonlinear/">Modeling Nonlinear Systems</a></li><li><a class="tocitem" href="../optimization/">Modeling Optimization Problems</a></li><li><a class="tocitem" href="../stochastic_diffeq/">Modeling with Stochasticity</a></li><li><a class="tocitem" href="../nonlinear_optimal_control/">Nonlinear Optimal Control</a></li></ul></li><li><span class="tocitem">ModelingToolkitize Tutorials</span><ul><li><a class="tocitem" href="../../mtkitize_tutorials/modelingtoolkitize/">Automatically Accelerating ODEProblem Code</a></li><li><a class="tocitem" href="../../mtkitize_tutorials/modelingtoolkitize_index_reduction/">Automated Index Reduction of DAEs</a></li></ul></li><li><span class="tocitem">Basics</span><ul><li><a class="tocitem" href="../../basics/AbstractSystem/">The AbstractSystem Interface</a></li><li><a class="tocitem" href="../../basics/ContextualVariables/">Contextual Variable Types</a></li><li><a class="tocitem" href="../../basics/Composition/">Composing Models and Building Reusable Components</a></li><li><a class="tocitem" href="../../basics/Validation/">Model Validation and Units</a></li><li><a class="tocitem" href="../../basics/DependencyGraphs/">Dependency Graphs</a></li><li><a class="tocitem" href="../../basics/FAQ/">Frequently Asked Questions</a></li></ul></li><li><span class="tocitem">System Types</span><ul><li><a class="tocitem" href="../../systems/ODESystem/">ODESystem</a></li><li><a class="tocitem" href="../../systems/SDESystem/">SDESystem</a></li><li><a class="tocitem" href="../../systems/JumpSystem/">JumpSystem</a></li><li><a class="tocitem" href="../../systems/NonlinearSystem/">NonlinearSystem</a></li><li><a class="tocitem" href="../../systems/OptimizationSystem/">OptimizationSystem</a></li><li><a class="tocitem" href="../../systems/ControlSystem/">ControlSystem</a></li><li><a class="tocitem" href="../../systems/ReactionSystem/">ReactionSystem</a></li><li><a class="tocitem" href="../../systems/PDESystem/">PDESystem</a></li></ul></li><li><a class="tocitem" href="../../comparison/">Comparison of ModelingToolkit vs Equation-Based Modeling Languages</a></li><li><a class="tocitem" href="../../internals/">Internal Details</a></li></ul><div class="docs-version-selector field has-addons"><div class="control"><span class="docs-label button is-static is-size-7">Version</span></div><div class="docs-selector control is-expanded"><div class="select is-fullwidth is-size-7"><select id="documenter-version-selector"></select></div></div></div></nav><div class="docs-main"><header class="docs-navbar"><nav class="breadcrumb"><ul class="is-hidden-mobile"><li><a class="is-disabled">Symbolic Modeling Tutorials</a></li><li class="is-active"><a href>Acausal Component-Based Modeling the RC Circuit</a></li></ul><ul class="is-hidden-tablet"><li class="is-active"><a href>Acausal Component-Based Modeling the RC Circuit</a></li></ul></nav><div class="docs-right"><a class="docs-edit-link" href="https://github.com/SciML/ModelingToolkit.jl/blob/master/docs/src/tutorials/acausal_components.md" title="Edit on GitHub"><span class="docs-icon fab"></span><span class="docs-label is-hidden-touch">Edit on GitHub</span></a><a class="docs-settings-button fas fa-cog" id="documenter-settings-button" href="#" title="Settings"></a><a class="docs-sidebar-button fa fa-bars is-hidden-desktop" id="documenter-sidebar-button" href="#"></a></div></header><article class="content" id="documenter-page"><h1 id="Acausal-Component-Based-Modeling-the-RC-Circuit"><a class="docs-heading-anchor" href="#Acausal-Component-Based-Modeling-the-RC-Circuit">Acausal Component-Based Modeling the RC Circuit</a><a id="Acausal-Component-Based-Modeling-the-RC-Circuit-1"></a><a class="docs-heading-anchor-permalink" href="#Acausal-Component-Based-Modeling-the-RC-Circuit" title="Permalink"></a></h1><p>In this tutorial we will build a hierarchical acausal component-based model of the RC circuit. The RC circuit is a simple example where we connect a resistor and a capacitor. <a href="https://en.wikipedia.org/wiki/Kirchhoff%27s_circuit_laws">Kirchoff's laws</a> are then applied to state equalities between currents and voltages. This specifies a differential-algebraic equation (DAE) system, where the algebraic equations are given by the constraints and equalities between different component variables. We then simplify this to an ODE by eliminating the equalities before solving. Let's see this in action.</p><h2 id="Copy-Paste-Example"><a class="docs-heading-anchor" href="#Copy-Paste-Example">Copy-Paste Example</a><a id="Copy-Paste-Example-1"></a><a class="docs-heading-anchor-permalink" href="#Copy-Paste-Example" title="Permalink"></a></h2><pre><code class="language-julia">using ModelingToolkit, Plots, DifferentialEquations
@parameters t
@connector function Pin(;name)
sts = @variables v(t)=1.0 i(t)=1.0
ODESystem(Equation[], t, sts, []; name=name)
end
function ModelingToolkit.connect(::Type{Pin}, ps...)
eqs = [
0 ~ sum(p->p.i, ps) # KCL
]
# KVL
for i in 1:length(ps)-1
push!(eqs, ps[i].v ~ ps[i+1].v)
end
return eqs
end
function Ground(;name)
@named g = Pin()
eqs = [g.v ~ 0]
compose(ODESystem(eqs, t, [], []; name=name), g)
end
function OnePort(;name)
@named p = Pin()
@named n = Pin()
sts = @variables v(t)=1.0 i(t)=1.0
eqs = [
v ~ p.v - n.v
0 ~ p.i + n.i
i ~ p.i
]
compose(ODESystem(eqs, t, sts, []; name=name), p, n)
end
function Resistor(;name, R = 1.0)
@named oneport = OnePort()
@unpack v, i = oneport
ps = @parameters R=R
eqs = [
v ~ i * R
]
extend(ODESystem(eqs, t, [], ps; name=name), oneport)
end
function Capacitor(;name, C = 1.0)
@named oneport = OnePort()
@unpack v, i = oneport
ps = @parameters C=C
D = Differential(t)
eqs = [
D(v) ~ i / C
]
extend(ODESystem(eqs, t, [], ps; name=name), oneport)
end
function ConstantVoltage(;name, V = 1.0)
@named oneport = OnePort()
@unpack v = oneport
ps = @parameters V=V
eqs = [
V ~ v
]
extend(ODESystem(eqs, t, [], ps; name=name), oneport)
end
R = 1.0
C = 1.0
V = 1.0
@named resistor = Resistor(R=R)
@named capacitor = Capacitor(C=C)
@named source = ConstantVoltage(V=V)
@named ground = Ground()
rc_eqs = [
connect(source.p, resistor.p)
connect(resistor.n, capacitor.p)
connect(capacitor.n, source.n, ground.g)
]
@named rc_model = compose(ODESystem(rc_eqs, t),
[resistor, capacitor, source, ground])
sys = structural_simplify(rc_model)
u0 = [
capacitor.v => 0.0
]
prob = ODAEProblem(sys, u0, (0, 10.0))
sol = solve(prob, Tsit5())
plot(sol)</code></pre><p><img src="https://user-images.githubusercontent.com/1814174/109416294-55184100-798b-11eb-9f05-766a793f0ba2.png" alt/></p><h2 id="Explanation"><a class="docs-heading-anchor" href="#Explanation">Explanation</a><a id="Explanation-1"></a><a class="docs-heading-anchor-permalink" href="#Explanation" title="Permalink"></a></h2><h3 id="Building-the-Component-Library"><a class="docs-heading-anchor" href="#Building-the-Component-Library">Building the Component Library</a><a id="Building-the-Component-Library-1"></a><a class="docs-heading-anchor-permalink" href="#Building-the-Component-Library" title="Permalink"></a></h3><p>For each of our components we use a Julia function which emits an <code>ODESystem</code>. At the top we start with defining the fundamental qualities of an electrical circuit component. At every input and output pin a circuit component has two values: the current at the pin and the voltage. Thus we define the <code>Pin</code> component (connector) to simply be the values there:</p><pre><code class="language-julia">@connector function Pin(;name)
sts = @variables v(t)=1.0 i(t)=1.0
ODESystem(Equation[], t, sts, []; name=name)
end</code></pre><p>Note that this is an incompletely specified ODESystem: it cannot be simulated on its own because the equations for <code>v(t)</code> and <code>i(t)</code> are unknown. Instead this just gives a common syntax for receiving this pair with some default values. Notice that in a component we define the <code>name</code> as a keyword argument: this is because later we will generate different <code>Pin</code> objects with different names to correspond to duplicates of this topology with unique variables. One can then construct a <code>Pin</code> like:</p><pre><code class="language-julia">Pin(name=:mypin1)</code></pre><p>or equivalently using the <code>@named</code> helper macro:</p><pre><code class="language-julia">@named mypin1 = Pin()</code></pre><p>Next we build our ground node. A ground node is just a pin that is connected to a constant voltage reservoir, typically taken to be <code>V=0</code>. Thus to define this component, we generate an <code>ODESystem</code> with a <code>Pin</code> subcomponent and specify that the voltage in such a <code>Pin</code> is equal to zero. This gives:</p><pre><code class="language-julia">function Ground(;name)
@named g = Pin()
eqs = [g.v ~ 0]
compose(ODESystem(eqs, t, [], []; name=name), g)
end</code></pre><p>Next we build a <code>OnePort</code>: an abstraction for all simple electrical component with two pins. The voltage difference between the positive pin and the negative pin is the voltage of the component, the current between two pins must sum to zero, and the current of the component equals to the current of the positive pin.</p><pre><code class="language-julia">function OnePort(;name)
@named p = Pin()
@named n = Pin()
sts = @variables v(t)=1.0 i(t)=1.0
eqs = [
v ~ p.v - n.v
0 ~ p.i + n.i
i ~ p.i
]
compose(ODESystem(eqs, t, sts, []; name=name), p, n)
end</code></pre><p>Next we build a resistor. A resistor is an object that has two <code>Pin</code>s, the positive and the negative pins, and follows Ohm's law: <code>v = i*r</code>. The voltage of the resistor is given as the voltage difference across the two pins while by conservation of charge we know that the current in must equal the current out, which means (no matter the direction of the current flow) the sum of the currents must be zero. This leads to our resistor equations:</p><pre><code class="language-julia">function Resistor(;name, R = 1.0)
@named oneport = OnePort()
@unpack v, i = oneport
ps = @parameters R=R
eqs = [
v ~ i * R
]
extend(ODESystem(eqs, t, [], ps; name=name, oneport)
end</code></pre><p>Notice that we have created this system with a default parameter <code>R</code> for the resistor's resistance. By doing so, if the resistance of this resistor is not overridden by a higher level default or overridden at <code>ODEProblem</code> construction time, this will be the value of the resistance. Also, note the use of <code>@unpack</code> and <code>extend</code>. For the <code>Resistor</code>, we want to simply inherit <code>OnePort</code>'s equations and states and extend them with a new equation. ModelingToolkit makes a new namespaced variable <code>oneport₊v(t)</code> when using the syntax <code>oneport.v</code>, and we can use <code>@unpack</code> avoid the namespacing.</p><p>Using our knowledge of circuits we similarly construct the <code>Capacitor</code>:</p><pre><code class="language-julia">function Capacitor(;name, C = 1.0)
@named oneport = OnePort()
@unpack v, i = oneport
ps = @parameters C=C
D = Differential(t)
eqs = [
D(v) ~ i / C
]
extend(ODESystem(eqs, t, [], ps; name=name), oneport)
end</code></pre><p>Now we want to build a constant voltage electrical source term. We can think of this as similarly being a two pin object, where the object itself is kept at a constant voltage, essentially generating the electrical current. We would then model this as:</p><pre><code class="language-julia">function ConstantVoltage(;name, V = 1.0)
@named oneport = OnePort()
@unpack v = oneport
ps = @parameters V=V
eqs = [
V ~ v
]
extend(ODESystem(eqs, t, [], ps; name=name), oneport)
end</code></pre><h3 id="Connecting-and-Simulating-Our-Electrical-Circuit"><a class="docs-heading-anchor" href="#Connecting-and-Simulating-Our-Electrical-Circuit">Connecting and Simulating Our Electrical Circuit</a><a id="Connecting-and-Simulating-Our-Electrical-Circuit-1"></a><a class="docs-heading-anchor-permalink" href="#Connecting-and-Simulating-Our-Electrical-Circuit" title="Permalink"></a></h3><p>Now we are ready to simulate our circuit. Let's build our four components: a <code>resistor</code>, <code>capacitor</code>, <code>source</code>, and <code>ground</code> term. For simplicity we will make all of our parameter values 1. This is done by:</p><pre><code class="language-julia">R = 1.0
C = 1.0
V = 1.0
@named resistor = Resistor(R=R)
@named capacitor = Capacitor(C=C)
@named source = ConstantVoltage(V=V)
@named ground = Ground()</code></pre><p>Next we have to define how we connect the circuit. Whenever two <code>Pin</code>s in a circuit are connected together, the system satisfies <a href="https://en.wikipedia.org/wiki/Kirchhoff%27s_circuit_laws">Kirchoff's laws</a>, i.e. that currents sum to zero and voltages across the pins are equal. Thus we will build a helper function <code>connect_pins</code> which implements these rules:</p><pre><code class="language-julia">function ModelingToolkit.connect(::Type{Pin}, ps...)
eqs = [
0 ~ sum(p->p.i, ps) # KCL
]
# KVL
for i in 1:length(ps)-1
push!(eqs, ps[i].v ~ ps[i+1].v)
end
return eqs
end</code></pre><p>Finally we will connect the pieces of our circuit together. Let's connect the positive pin of the resistor to the source, the negative pin of the resistor to the capacitor, and the negative pin of the capacitor to a junction between the source and the ground. This would mean our connection equations are:</p><pre><code class="language-julia">rc_eqs = [
connect(source.p, resistor.p)
connect(resistor.n, capacitor.p)
connect(capacitor.n, source.n, ground.g)
]</code></pre><p>Finally we build our four component model with these connection rules:</p><pre><code class="language-julia">@named rc_model = compose(ODESystem(rc_eqs, t)
[resistor, capacitor, source, ground])</code></pre><p>Note that we can also specify the subsystems in a vector. This model is acasual because we have not specified anything about the causality of the model. We have simply specified what is true about each of the variables. This forms a system of differential-algebraic equations (DAEs) which define the evolution of each state of the system. The equations are:</p><pre><code class="language-julia">equations(rc_model)
20-element Vector{Equation}:
0 ~ resistor₊p₊i(t) + source₊p₊i(t)
source₊p₊v(t) ~ resistor₊p₊v(t)
0 ~ capacitor₊p₊i(t) + resistor₊n₊i(t)
resistor₊n₊v(t) ~ capacitor₊p₊v(t)
0 ~ capacitor₊n₊i(t) + ground₊g₊i(t) + source₊n₊i(t)
capacitor₊n₊v(t) ~ source₊n₊v(t)
source₊n₊v(t) ~ ground₊g₊v(t)
resistor₊v(t) ~ resistor₊p₊v(t) - resistor₊n₊v(t)
0 ~ resistor₊n₊i(t) + resistor₊p₊i(t)
resistor₊i(t) ~ resistor₊p₊i(t)
resistor₊v(t) ~ resistor₊R*resistor₊i(t)
capacitor₊v(t) ~ capacitor₊p₊v(t) - capacitor₊n₊v(t)
0 ~ capacitor₊n₊i(t) + capacitor₊p₊i(t)
capacitor₊i(t) ~ capacitor₊p₊i(t)
Differential(t)(capacitor₊v(t)) ~ capacitor₊i(t)*(capacitor₊C^-1)
source₊v(t) ~ source₊p₊v(t) - source₊n₊v(t)
0 ~ source₊n₊i(t) + source₊p₊i(t)
source₊i(t) ~ source₊p₊i(t)
source₊V ~ source₊v(t)
ground₊g₊v(t) ~ 0</code></pre><p>the states are:</p><pre><code class="language-julia">states(rc_model)
20-element Vector{Term{Real, Base.ImmutableDict{DataType, Any}}}:
source₊p₊i(t)
resistor₊p₊i(t)
source₊p₊v(t)
resistor₊p₊v(t)
capacitor₊p₊i(t)
resistor₊n₊i(t)
resistor₊n₊v(t)
capacitor₊p₊v(t)
source₊n₊i(t)
capacitor₊n₊i(t)
ground₊g₊i(t)
capacitor₊n₊v(t)
source₊n₊v(t)
ground₊g₊v(t)
resistor₊v(t)
resistor₊i(t)
capacitor₊v(t)
capacitor₊i(t)
source₊v(t)
source₊i(t)</code></pre><p>and the parameters are:</p><pre><code class="language-julia">parameters(rc_model)
3-element Vector{Any}:
resistor₊R
capacitor₊C
source₊V</code></pre><h2 id="Simplifying-and-Solving-this-System"><a class="docs-heading-anchor" href="#Simplifying-and-Solving-this-System">Simplifying and Solving this System</a><a id="Simplifying-and-Solving-this-System-1"></a><a class="docs-heading-anchor-permalink" href="#Simplifying-and-Solving-this-System" title="Permalink"></a></h2><p>This system could be solved directly as a DAE using <a href="https://diffeq.sciml.ai/stable/solvers/dae_solve/">one of the DAE solvers from DifferentialEquations.jl</a>. However, let's take a second to symbolically simplify the system before doing the solve. Although we can use ODE solvers that handles mass matrices to solve the above system directly, we want to run the <code>structural_simplify</code> function first, as it eliminates many unnecessary variables to build the leanest numerical representation of the system. Let's see what it does here:</p><pre><code class="language-julia">sys = structural_simplify(rc_model)
equations(sys)
2-element Vector{Equation}:
0 ~ capacitor₊v(t) + resistor₊R*resistor₊i(t) - source₊V
Differential(t)(capacitor₊v(t)) ~ resistor₊i(t)*(capacitor₊C^-1)</code></pre><pre><code class="language-julia">states(sys)
2-element Vector{Any}:
capacitor₊v(t)
capacitor₊p₊i(t)</code></pre><p>After structural simplification we are left with a system of only two equations with two state variables. One of the equations is a differential equation while the other is an algebraic equation. We can then give the values for the initial conditions of our states and solve the system by converting it to an ODEProblem in mass matrix form and solving it with an <a href="https://diffeq.sciml.ai/stable/solvers/dae_solve/#OrdinaryDiffEq.jl-(Mass-Matrix)">ODEProblem mass matrix DAE solver</a>. This is done as follows:</p><pre><code class="language-julia">u0 = [
capacitor.v => 0.0
capacitor.p.i => 0.0
]
prob = ODEProblem(sys, u0, (0, 10.0))
sol = solve(prob, Rodas4())
plot(sol)</code></pre><p><img src="https://user-images.githubusercontent.com/1814174/109416295-55184100-798b-11eb-96d1-5bb7e40135ba.png" alt/></p><p>Since we have run <code>structural_simplify</code>, MTK can numerically solve all the unreduced algebraic equations numerically using the <code>ODAEProblem</code> (note the letter <code>A</code>):</p><pre><code class="language-julia">u0 = [
capacitor.v => 0.0
]
prob = ODAEProblem(sys, u0, (0, 10.0))
sol = solve(prob, Rodas4())
plot(sol)</code></pre><p><img src="https://user-images.githubusercontent.com/1814174/109416294-55184100-798b-11eb-9f05-766a793f0ba2.png" alt/></p><p>Notice that this solves the whole system by only solving for one variable!</p><p>However, what if we wanted to plot the timeseries of a different variable? Do not worry, that information was not thrown away! Instead, transformations like <code>structural_simplify</code> simply change state variables into <code>observed</code> variables. Let's see what our observed variables are:</p><pre><code class="language-julia">observed(sys)
18-element Vector{Equation}:
capacitor₊i(t) ~ resistor₊i(t)
ground₊g₊i(t) ~ 0.0
source₊n₊i(t) ~ resistor₊i(t)
source₊i(t) ~ -resistor₊i(t)
source₊p₊i(t) ~ -resistor₊i(t)
capacitor₊n₊i(t) ~ -resistor₊i(t)
resistor₊n₊v(t) ~ capacitor₊v(t)
resistor₊n₊i(t) ~ -resistor₊i(t)
resistor₊p₊i(t) ~ resistor₊i(t)
capacitor₊p₊i(t) ~ resistor₊i(t)
capacitor₊p₊v(t) ~ capacitor₊v(t)
capacitor₊n₊v(t) ~ 0.0
source₊n₊v(t) ~ 0.0
ground₊g₊v(t) ~ 0.0
source₊v(t) ~ source₊V
source₊p₊v(t) ~ source₊v(t)
resistor₊p₊v(t) ~ source₊v(t)
resistor₊v(t) ~ source₊v(t) - capacitor₊v(t)</code></pre><p>These are explicit algebraic equations which can then be used to reconstruct the required variables on the fly. This leads to dramatic computational savings because implicitly solving an ODE scales like O(n^3), so making there be as few states as possible is good!</p><p>The solution object can be accessed via its symbols. For example, let's retrieve the voltage of the resistor over time:</p><pre><code class="language-julia">sol[resistor.v]</code></pre><p>or we can plot the timeseries of the resistor's voltage:</p><pre><code class="language-julia">plot(sol, vars=[resistor.v])</code></pre></article><nav class="docs-footer"><a class="docs-footer-prevpage" href="../ode_modeling/">« Composing Ordinary Differential Equations</a><a class="docs-footer-nextpage" href="../higher_order/">Automatic Transformation of Nth Order ODEs to 1st Order ODEs »</a><div class="flexbox-break"></div><p class="footer-message">Powered by <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> and the <a href="https://julialang.org/">Julia Programming Language</a>.</p></nav></div><div class="modal" id="documenter-settings"><div class="modal-background"></div><div class="modal-card"><header class="modal-card-head"><p class="modal-card-title">Settings</p><button class="delete"></button></header><section class="modal-card-body"><p><label class="label">Theme</label><div class="select"><select id="documenter-themepicker"><option value="documenter-light">documenter-light</option><option value="documenter-dark">documenter-dark</option></select></div></p><hr/><p>This document was generated with <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> version 0.27.3 on <span class="colophon-date" title="Tuesday 13 July 2021 11:09">Tuesday 13 July 2021</span>. Using Julia version 1.6.1.</p></section><footer class="modal-card-foot"></footer></div></div></div></body></html>