|
|
| Line 84: |
Line 84: |
| | plot(subs(soln2, g = 9.81, y0 = 5, v0 = 10, y(t)), t = 0 .. 2) | | plot(subs(soln2, g = 9.81, y0 = 5, v0 = 10, y(t)), t = 0 .. 2) |
| | </syntaxhighlight> | | </syntaxhighlight> |
| − | and can of course add other plotting options as needed. Note that ''soln2'' is a single equation, not a collection. If you end up getting a collection of results, you may need to de-bracket them. | + | and can of course add other plotting options as needed. Once again, note that ''soln2'' is a single equation, not a collection, so there is no need to de-bracket. |
| − | | |
| − | | |
| − | <!--
| |
| − | == Basic Symbolic Example ==
| |
| − | For some simple systems of differential equations, SymPy can actually solve with symbolic coefficients and undefined forcing functions. The result will be an unevaluated integral, but still - pretty impressive! Here is an example that solves $$a\,\frac{dy(t)}{dt} + b\,y(t)=x(t)$$ assuming $$y(t_0)=y_0$$. If you want to follow along, the code is in the <code>dsolve_2</code> notebook in the [https://drive.google.com/drive/folders/10Y1kzpu-DTTRwm1mYtdgTwpXF2zKpmoL?usp=sharing Basic Dsolve] Google Drive folder.
| |
| − | <syntaxhighlight lang=python>
| |
| − | import sympy as sym
| |
| − | sym.init_printing(use_latex=True)
| |
| − | a, b, t, y0, t0 = sym.var('a, b, t, y0, t0')
| |
| − | x = sym.Function('x')
| |
| − | y = sym.Function('y')
| |
| − | eqn1 = sym.Eq(a*sym.Derivative(y(t), t) + b * y(t), x(t))
| |
| − | eqns = [eqn1]
| |
| − | soln = sym.dsolve(eqns, [y(t)], ics = {y(t0): y0})
| |
| − | soln
| |
| − | </syntaxhighlight>
| |
| − | The solution is:
| |
| − | <center>
| |
| − | $$
| |
| − | \left[ y{\left(t \right)} = \left(y_{0} e^{\frac{b t_{0}}{a}} - \frac{\int x{\left(t_{0} \right)} e^{\frac{b t_{0}}{a}}\, dt_{0}}{a}\right) e^{- \frac{b t}{a}} + \frac{e^{- \frac{b t}{a}} \int x{\left(t \right)} e^{\frac{b t}{a}}\, dt}{a}\right]
| |
| − | $$
| |
| − | </center>
| |
| − | === Making Substitutions ===
| |
| − | Making substitutions into lists of equations is a little complicated - you need to go into each element in the list and make the substitutions there. Given that, there is a helper function in <code>sym_helper</code> called <code>makesubseqnlist(EQN LIST, VARS TO SUB, VALS TO SUB)</code> that will take a list of equations, a list of variables (or functions), and a list of values (or functions). The code for <code>makesubseqnlist(EQN LIST, VARS TO SUB, VALS TO SUB)</code> is
| |
| − | | |
| − | <syntaxhighlight lang=python>
| |
| − | def makesubseqnlist(equations, vartosub, valtosub):
| |
| − | sublist = list(zip(vartosub, valtosub))
| |
| − | esub = []
| |
| − | for eqn in equations:
| |
| − | esub.append(eqn.subs(sublist))
| |
| − | return esub
| |
| − | </syntaxhighlight>
| |
| − | | |
| − | and you can get by by adding the following to near the top of your notebook or script:
| |
| − | <syntaxhighlight lang=python>
| |
| − | try:
| |
| − | from sym_helper import *
| |
| − | except: # Get sym_helper.py from GitHub
| |
| − | import requests
| |
| − | url = 'https://tinyurl.com/mrgdemos/SymPy_Files/sym_helper.py'
| |
| − | open('sym_helper.py', 'w').write(requests.get(url).text)
| |
| − | from sym_helper import *
| |
| − | </syntaxhighlight>
| |
| − | | |
| − | If we specifically want to set $$a=2$$, $$b=3$$, and $$x(t)=4\,\cos(5t)$$, and then solve that where $$y(0)=6$$, we can:
| |
| − | <syntaxhighlight lang=python>
| |
| − | eqns2 = makesubseqnlist(eqns, [a, b, x(t)], [2, 3, 4*sym.cos(5*t)])
| |
| − | soln2 = sym.dsolve(eqns2, [y(t)], ics = {y(0): 6})
| |
| − | </syntaxhighlight>
| |
| − | yields:
| |
| − | <center>
| |
| − | $$
| |
| − | \left[ y{\left(t \right)} = \frac{40 \sin{\left(5 t \right)}}{109} + \frac{12 \cos{\left(5 t \right)}}{109} + \frac{642 e^{- \frac{3 t}{2}}}{109}\right]
| |
| − | $$
| |
| − | </center>
| |
| − | | |
| − | == Making a Dictionary ==
| |
| − | For linear algebra, the results came in a dictionary that had each of the unknowns as keys and their solutions as variables. Since <code>dsolve</code> produces a list, versus a ductionary, it will be useful to convert that list to a dictionary for making substitutions later. <code>sym_helper</code> has a <code>make_dict(LIST)</code> command that will turn a list of equations into a dictionary. The code for <code>make_dict</code> is:
| |
| − | <syntaxhighlight lang=python>
| |
| − | def makedict(sollist):
| |
| − | d = {}
| |
| − | for eqn in sollist:
| |
| − | d[eqn.lhs] = eqn.rhs
| |
| − | return d
| |
| − | </syntaxhighlight>
| |
| − | and so the code you can add to your script is:
| |
| − | <syntaxhighlight lang=python>
| |
| − | def makedict(sollist):
| |
| − | soln = makedict(solnlist)
| |
| − | </syntaxhighlight>
| |
| − | and, as long as <code>sym_helper</code> has been imported, you can run it by writing code like:
| |
| − | <syntaxhighlight lang=python>
| |
| − | solndict = makedict(soln2)
| |
| − | </syntaxhighlight>
| |
| − | which yields:
| |
| − | <center>
| |
| − | $$
| |
| − | \left\{ y{\left(t \right)} : \frac{40 \sin{\left(5 t \right)}}{109} + \frac{12 \cos{\left(5 t \right)}}{109} + \frac{642 e^{- \frac{3 t}{2}}}{109}\right\}
| |
| − | | |
| − | $$
| |
| − | </center>
| |
| − | | |
| − | == Plotting Results ==
| |
| − | Plotting will require using <code>lambdify</code> to turn one or more functions into a lambda function and then setting up your plot. Here is the code for this example that does all those things (making a plot that runs from $$t$$=0 to 10 sec); note the additional required imports:
| |
| − | <syntaxhighlight lang=python>
| |
| − | import matplotlib.pyplot as plt
| |
| − | import numpy as np
| |
| − | yfun = sym.lambdify(t, solndict[y(t)])
| |
| − | fig = plt.figure(num=1, clear=True)
| |
| − | ax = fig.add_subplot(1, 1, 1)
| |
| − | tn = np.linspace(0, 10, 1000)
| |
| − | ax.plot(tn, yfun(tn))
| |
| − | </syntaxhighlight>
| |
| − | | |
| − | == Lessons Learned ==
| |
| − | * dsolve is good for up to two coupled ODEs. Once there are three or more...not so much (even if all the coefficients are ints or floats). Once there are three or more (linear) coupled differential equations, consider going the Laplace route.
| |
| − | * If your equations are a combination of differential equations a non-differential equations, dsolve expects all the differential equations to come first.
| |
| − | * If you are using certain functions from SymPy that aren't exactly the same in NumPy, when you lambdify things you will need to provide the translation
| |
| − | == Examples ==
| |
| − | * See [[SymPy/Differential Equations/RC Example]] for a complete example; this will reference IVP_RCDemo in the Google Drive folder.
| |
| − | * See [[SymPy/CAD]] for a list of SymPy examples from a free textbook. The ones most relevant to this page are from Chapter 5.
| |
| − | * There are several other files in the [https://drive.google.com/drive/folders/10Y1kzpu-DTTRwm1mYtdgTwpXF2zKpmoL?usp=sharing Basic Dsolve] Google Drive folder:
| |
| − | | |
| − | ** Coupled ODEs (just to see how SymPy handles those!)
| |
| − | *** Coupled2
| |
| − | *** Coupled3
| |
| − | *** deq_coupled: symbolically solve two coupled first-order ODEs; substitute numbers later to make graph. Needs some editing - does not completely follow the general template.
| |
| − | *** deq_coupled_num: symbolically solve two coupled first-order ODEs with numbers substituted in before solving. Needs some editing - does not completely follow the general template.
| |
| − | ** First-order differential equations:
| |
| − | *** IVP_demo1st: symbolic solution to 1st order ODE with a constant forcing function; substitutes numbers later to make graph
| |
| − | *** IVP_demo1stfun: symbolic solution to 1st order ODE with a non-constant forcing function; substitutes numbers later to make graph (including forcing function on graph)
| |
| − | *** IVP_demo1st_pulse: symbolic solution to 1st order ODE with a non-constant forcing function comprised of step and ramp functions; substitutes numbers later to make graph (including forcing function on graph)
| |
| − | **** This one shows extra code needed to translate step and ramp functions into a lambdify function
| |
| − | ** Second-order differential equations:
| |
| − | *** IVP_demo2ndcon: symbolic solution to 2nd order ODE with a constant forcing function; substitutes numbers later to make graph. Needs some editing - does not completely follow the general template.
| |
| − | *** IVP_demo2ndfun: symbolic solution to 2nd order ODE with a non-constant forcing function; substitutes numbers later to make graph. Needs some editing - does not completely follow the general template.
| |
| − | -->
| |
Introduction
This page focuses on using Maple to find both the symbolic and the numeric solutions to differential equations obtained from electric circuits.
Note: There is a...decidedly more complicated explanation of these things at Maple/Differential Equations/Old; the goal for Spring 2024 and beyond is to keep things simpler.
Very Basic Example
Imagine you have the equation $$2\frac{dy(t)}{dt} + 3 y(t) = 4$$ with the initial condition $$y(0)=5$$, and you want to solve for $$y(t)$$. You can do this with Maple as follows:
Initialization
As always, it is a good idea to include the restart command in your worksheet:
Define Equation
In the same way that you were able to assign a linear algebra expression to a variable, you can do the same with a differential equation. The key is to note that the Maple diff command can be used to calculate or represent a derivative. You will need to explicitly let Maple know that your variable is a function (in our case, a function of $$t$$) by including that parameter with the variable. Given that, you can store the differential equation in a variable called eqn with:
eqn1:=2*diff(y(t), t)+3*y(t)=4
Solve Equation
Solving a system of differential equations is also similar to solving a system of linear algebra equations - the main differences are that you will use dsolve instead of solve, you must continue to use $$y(t)$$ instead of just $$y$$, and you may end up needing to add some initial conditions. The code
soln1:=dsolve([eqn1], [y(t)])
will result in a solution of:
\(\mathit{soln1}:=\left\{y\! \left(t\right)=\frac{4}{3}+{\mathrm e}^{-\frac{3 t}{2}} \textit{_}\mathit{C1}\right\}\)
Initial Condition
To incorporate initial conditions, you will give the dsolve command information about the value of the variable (or, for higher order differential equations, the value and values of the derivatives of the variable). For example, to solve our sample equation with $$y(0)=5$$, you will include the initial condition by adding y(0)=5 to the equations:
soln2:=dsolve([eqn1, y(0)=5], [y(t)])
will produce
\(
\mathit{soln2}:=y\! \left(t\right)=\frac{4}{3}+\frac{11 {\mathrm e}^{-\frac{3 t}{2}}}{3}
\)
Note that the initial condition does not have to be at time 0; if you know that $$y(6)=7$$, you can use that as well:
soln3:=dsolve([eqn1, y(6)=7], [y(t)])
will produce
\(
\mathit{soln3}:=y\! \left(t\right)=\frac{4}{3}+\frac{17 {\mathrm e}^{-\frac{3 t}{2}}}{3}
\)
Making a Plot
Once you have a solution or set of solutions, you can plot them using subs. For instance, to plot $$y(t)$$ in soln2 above (where we set $$y(0)=5$$) you can use:
plot(subs(soln2, y(t)), t = 0 .. 5)
and can of course add other plotting options as needed. Note that soln2 is a single equation, not a collection. If you end up getting a collection of results, you may need to de-bracket them.
Second Order Example
Now imagine that you want to solve for $$y(t)$$ in
\(
\frac{d^2y(t)}{dt^2}=-g\)
We can solve this using all symbols or numbers, and we can solve with or without initial conditions.
Initialization
Define Equation
eqn1 := diff(y(t), t, t) = -g
Solve Equation / Initial Conditions
Without initial conditions, you can use
soln1 := dsolve([eqn1], [y(t)])
to get:
\(\mathit{soln1}:=\left\{y\! \left(t\right)=-\frac{1}{2} g \,t^{2}+\textit{_}\mathit{C1} t+\textit{_}\mathit{C2}\right\}\)
You can also put in symbolic initial conditions. To put in a derivative condition, use the format D^n(var)(time)=val to establish a condition for the nth derivative of var at time time. For instance, if you know some initial velocity and position, we can use:
soln2 := dsolve([eqn1, y(0) = y0, D(y)(0) = v0], [y(t)])
to get:
\(\mathit{soln2}:=y\! \left(t\right)=-\frac{1}{2} g \,t^{2}+\mathit{v0} t+\mathit{y0}\)
Making a Plot
Since your answer has symbols, you will need to replace them with numbers before plotting:
plot(subs(soln2, g = 9.81, y0 = 5, v0 = 10, y(t)), t = 0 .. 2)
and can of course add other plotting options as needed. Once again, note that soln2 is a single equation, not a collection, so there is no need to de-bracket.
