In [1]:
%load_ext autoreload
%autoreload 2
import sys
import pylab as pl
sys.path.insert(0, '..')
import pymoca.parser
import pymoca.backends.sympy.generator as generator
import pylab as pl
import sympy
import control
sympy.init_printing()
%load_ext autoreload
%autoreload 2
%matplotlib inline
The autoreload extension is already loaded. To reload it, use:
%reload_ext autoreload
In [2]:
modelica_src = '''
model System
input Real u;
output Real x(start=1), v_x(start=1);
Spring spring;
Damper damper;
equation
spring.x = x;
damper.v = v_x;
der(x) = v_x;
der(v_x) = spring.f + damper.f - u;
end System;
model Spring
Real x "displacement";
Real f "force";
parameter Real k = 2.0 "spring constant";
equation
f = -k*x;
end Spring;
model Damper
Real v "velocity";
Real f "force";
parameter Real c = 0.2 "damping constant";
equation
f = -c*v;
end Damper;
'''
In [30]:
ast = pymoca.parser.parse(modelica_src)
ast.classes['System'].symbols.keys()
Out[30]:
odict_keys(['u', 'x', 'v_x', 'spring', 'damper'])
In [32]:
ast
Out[32]:
{
"annotation": [],
"classes": {
"Damper": {
"annotation": [],
"classes": {},
"comment": "",
"encapsulated": false,
"equations": [
{
"comment": "",
"left": {
"child": [],
"indices": [],
"name": "f"
},
"right": {
"operands": [
{
"operands": [
{
"child": [],
"indices": [],
"name": "c"
}
],
"operator": "-"
},
{
"child": [],
"indices": [],
"name": "v"
}
],
"operator": "*"
}
}
],
"extends": [],
"final": false,
"functions": {},
"imports": [],
"initial_equations": [],
"initial_statements": [],
"name": "Damper",
"partial": false,
"statements": [],
"symbols": {
"c": {
"class_modification": {
"arguments": [
{
"redeclare": false,
"value": {
"component": {
"child": [],
"indices": [],
"name": "value"
},
"modifications": [
{
"value": 0.2
}
]
}
}
]
},
"comment": "damping constant",
"dimensions": [
{
"value": 1
}
],
"final": false,
"fixed": {
"value": false
},
"id": 0,
"inner": false,
"max": {
"value": null
},
"min": {
"value": null
},
"name": "c",
"nominal": {
"value": null
},
"order": 10,
"outer": false,
"prefixes": [
"parameter"
],
"redeclare": false,
"start": {
"value": 0
},
"type": {
"child": [],
"indices": [],
"name": "Real"
},
"value": {
"value": null
},
"visibility": "private"
},
"f": {
"class_modification": null,
"comment": "force",
"dimensions": [
{
"value": 1
}
],
"final": false,
"fixed": {
"value": false
},
"id": 0,
"inner": false,
"max": {
"value": null
},
"min": {
"value": null
},
"name": "f",
"nominal": {
"value": null
},
"order": 9,
"outer": false,
"prefixes": [],
"redeclare": false,
"start": {
"value": 0
},
"type": {
"child": [],
"indices": [],
"name": "Real"
},
"value": {
"value": null
},
"visibility": "private"
},
"v": {
"class_modification": null,
"comment": "velocity",
"dimensions": [
{
"value": 1
}
],
"final": false,
"fixed": {
"value": false
},
"id": 0,
"inner": false,
"max": {
"value": null
},
"min": {
"value": null
},
"name": "v",
"nominal": {
"value": null
},
"order": 8,
"outer": false,
"prefixes": [],
"redeclare": false,
"start": {
"value": 0
},
"type": {
"child": [],
"indices": [],
"name": "Real"
},
"value": {
"value": null
},
"visibility": "private"
}
},
"type": "model"
},
"Spring": {
"annotation": [],
"classes": {},
"comment": "",
"encapsulated": false,
"equations": [
{
"comment": "",
"left": {
"child": [],
"indices": [],
"name": "f"
},
"right": {
"operands": [
{
"operands": [
{
"child": [],
"indices": [],
"name": "k"
}
],
"operator": "-"
},
{
"child": [],
"indices": [],
"name": "x"
}
],
"operator": "*"
}
}
],
"extends": [],
"final": false,
"functions": {},
"imports": [],
"initial_equations": [],
"initial_statements": [],
"name": "Spring",
"partial": false,
"statements": [],
"symbols": {
"f": {
"class_modification": null,
"comment": "force",
"dimensions": [
{
"value": 1
}
],
"final": false,
"fixed": {
"value": false
},
"id": 0,
"inner": false,
"max": {
"value": null
},
"min": {
"value": null
},
"name": "f",
"nominal": {
"value": null
},
"order": 6,
"outer": false,
"prefixes": [],
"redeclare": false,
"start": {
"value": 0
},
"type": {
"child": [],
"indices": [],
"name": "Real"
},
"value": {
"value": null
},
"visibility": "private"
},
"k": {
"class_modification": {
"arguments": [
{
"redeclare": false,
"value": {
"component": {
"child": [],
"indices": [],
"name": "value"
},
"modifications": [
{
"value": 2.0
}
]
}
}
]
},
"comment": "spring constant",
"dimensions": [
{
"value": 1
}
],
"final": false,
"fixed": {
"value": false
},
"id": 0,
"inner": false,
"max": {
"value": null
},
"min": {
"value": null
},
"name": "k",
"nominal": {
"value": null
},
"order": 7,
"outer": false,
"prefixes": [
"parameter"
],
"redeclare": false,
"start": {
"value": 0
},
"type": {
"child": [],
"indices": [],
"name": "Real"
},
"value": {
"value": null
},
"visibility": "private"
},
"x": {
"class_modification": null,
"comment": "displacement",
"dimensions": [
{
"value": 1
}
],
"final": false,
"fixed": {
"value": false
},
"id": 0,
"inner": false,
"max": {
"value": null
},
"min": {
"value": null
},
"name": "x",
"nominal": {
"value": null
},
"order": 5,
"outer": false,
"prefixes": [],
"redeclare": false,
"start": {
"value": 0
},
"type": {
"child": [],
"indices": [],
"name": "Real"
},
"value": {
"value": null
},
"visibility": "private"
}
},
"type": "model"
},
"System": {
"annotation": [],
"classes": {},
"comment": "",
"encapsulated": false,
"equations": [
{
"comment": "",
"left": {
"child": [
{
"child": [],
"indices": [],
"name": "x"
}
],
"indices": [],
"name": "spring"
},
"right": {
"child": [],
"indices": [],
"name": "x"
}
},
{
"comment": "",
"left": {
"child": [
{
"child": [],
"indices": [],
"name": "v"
}
],
"indices": [],
"name": "damper"
},
"right": {
"child": [],
"indices": [],
"name": "v_x"
}
},
{
"comment": "",
"left": {
"operands": [
{
"child": [],
"indices": [],
"name": "x"
}
],
"operator": "der"
},
"right": {
"child": [],
"indices": [],
"name": "v_x"
}
},
{
"comment": "",
"left": {
"operands": [
{
"child": [],
"indices": [],
"name": "v_x"
}
],
"operator": "der"
},
"right": {
"operands": [
{
"operands": [
{
"child": [
{
"child": [],
"indices": [],
"name": "f"
}
],
"indices": [],
"name": "spring"
},
{
"child": [
{
"child": [],
"indices": [],
"name": "f"
}
],
"indices": [],
"name": "damper"
}
],
"operator": "+"
},
{
"child": [],
"indices": [],
"name": "u"
}
],
"operator": "-"
}
}
],
"extends": [],
"final": false,
"functions": {},
"imports": [],
"initial_equations": [],
"initial_statements": [],
"name": "System",
"partial": false,
"statements": [],
"symbols": {
"damper": {
"class_modification": null,
"comment": "",
"dimensions": [
{
"value": 1
}
],
"final": false,
"fixed": {
"value": false
},
"id": 0,
"inner": false,
"max": {
"value": null
},
"min": {
"value": null
},
"name": "damper",
"nominal": {
"value": null
},
"order": 4,
"outer": false,
"prefixes": [],
"redeclare": false,
"start": {
"value": 0
},
"type": {
"child": [],
"indices": [],
"name": "Damper"
},
"value": {
"value": null
},
"visibility": "private"
},
"spring": {
"class_modification": null,
"comment": "",
"dimensions": [
{
"value": 1
}
],
"final": false,
"fixed": {
"value": false
},
"id": 0,
"inner": false,
"max": {
"value": null
},
"min": {
"value": null
},
"name": "spring",
"nominal": {
"value": null
},
"order": 3,
"outer": false,
"prefixes": [],
"redeclare": false,
"start": {
"value": 0
},
"type": {
"child": [],
"indices": [],
"name": "Spring"
},
"value": {
"value": null
},
"visibility": "private"
},
"u": {
"class_modification": null,
"comment": "",
"dimensions": [
{
"value": 1
}
],
"final": false,
"fixed": {
"value": false
},
"id": 0,
"inner": false,
"max": {
"value": null
},
"min": {
"value": null
},
"name": "u",
"nominal": {
"value": null
},
"order": 0,
"outer": false,
"prefixes": [
"input"
],
"redeclare": false,
"start": {
"value": 0
},
"type": {
"child": [],
"indices": [],
"name": "Real"
},
"value": {
"value": null
},
"visibility": "private"
},
"v_x": {
"class_modification": {
"arguments": [
{
"redeclare": false,
"value": {
"component": {
"child": [],
"indices": [],
"name": "start"
},
"modifications": [
{
"value": 1.0
}
]
}
}
]
},
"comment": "",
"dimensions": [
{
"value": 1
}
],
"final": false,
"fixed": {
"value": false
},
"id": 0,
"inner": false,
"max": {
"value": null
},
"min": {
"value": null
},
"name": "v_x",
"nominal": {
"value": null
},
"order": 2,
"outer": false,
"prefixes": [
"output"
],
"redeclare": false,
"start": {
"value": 0
},
"type": {
"child": [],
"indices": [],
"name": "Real"
},
"value": {
"value": null
},
"visibility": "private"
},
"x": {
"class_modification": {
"arguments": [
{
"redeclare": false,
"value": {
"component": {
"child": [],
"indices": [],
"name": "start"
},
"modifications": [
{
"value": 1.0
}
]
}
}
]
},
"comment": "",
"dimensions": [
{
"value": 1
}
],
"final": false,
"fixed": {
"value": false
},
"id": 0,
"inner": false,
"max": {
"value": null
},
"min": {
"value": null
},
"name": "x",
"nominal": {
"value": null
},
"order": 1,
"outer": false,
"prefixes": [
"output"
],
"redeclare": false,
"start": {
"value": 0
},
"type": {
"child": [],
"indices": [],
"name": "Real"
},
"value": {
"value": null
},
"visibility": "private"
}
},
"type": "model"
}
},
"comment": "",
"encapsulated": false,
"equations": [],
"extends": [],
"final": false,
"functions": {},
"imports": [],
"initial_equations": [],
"initial_statements": [],
"name": null,
"partial": false,
"statements": [],
"symbols": {},
"type": ""
}
In [7]:
src_code = generator.generate(ast, 'System')
print(src_code)
# do not edit, generated by pymoca
from __future__ import print_function, division
import sympy
import sympy.physics.mechanics as mech
from pymoca.backends.sympy.runtime import OdeModel
from sympy import sin, cos, tan
class System(OdeModel):
def __init__(self):
super(System, self).__init__()
# states
x, v_x = mech.dynamicsymbols('x, v_x')
self.x = sympy.Matrix([x, v_x])
self.x0 = {
x : 1.0,
v_x : 1.0,
}
# variables
spring__x, spring__f, damper__v, damper__f = mech.dynamicsymbols('spring.x, spring.f, damper.v, damper.f')
self.v = sympy.Matrix([spring__x, spring__f, damper__v, damper__f])
# constants
self.c = sympy.Matrix([])
self.c0 = {
}
# parameters
spring__k, damper__c = sympy.symbols('spring.k, damper.c')
self.p = sympy.Matrix([spring__k, damper__c])
self.p0 = {
spring__k : 2.0,
damper__c : 0.2,
}
# inputs
u = mech.dynamicsymbols('u')
self.u = sympy.Matrix([u])
self.u0 = {
u : 0,
}
# outputs
x, v_x = mech.dynamicsymbols('x, v_x')
self.y = sympy.Matrix([x, v_x])
# equations
self.eqs = [
spring__f - (- spring__k * spring__x),
damper__f - (- damper__c * damper__v),
spring__x - (x),
damper__v - (v_x),
(x).diff(self.t) - (v_x),
(v_x).diff(self.t) - (spring__f + damper__f - u),
]
self.compute_fg()
In [8]:
exec(src_code)
model = System()
model
Out[8]:
{'t': t, 'x': Matrix([
[ x(t)],
[v_x(t)]]), 'u': Matrix([[u(t)]]), 'y': Matrix([
[ x(t)],
[v_x(t)]]), 'p': Matrix([
[spring.k],
[damper.c]]), 'c': Matrix(0, 0, []), 'v': Matrix([
[spring.x(t)],
[spring.f(t)],
[damper.v(t)],
[damper.f(t)]]), 'x0': {x(t): 1.0, v_x(t): 1.0}, 'u0': {u(t): 0}, 'p0': {spring.k: 2.0, damper.c: 0.2}, 'c0': {}, 'eqs': [spring.k*spring.x(t) + spring.f(t), damper.c*damper.v(t) + damper.f(t), spring.x(t) - x(t), damper.v(t) - v_x(t), -v_x(t) + Derivative(x(t), t), -damper.f(t) - spring.f(t) + u(t) + Derivative(v_x(t), t)], 'f': Matrix([
[ v_x(t)],
[-damper.c*v_x(t) - spring.k*x(t) - u(t)]]), 'g': Matrix([
[ x(t)],
[v_x(t)]])}
In [9]:
t = sympy.Symbol('t')
model.x.diff(t)
Out[9]:
$$\left[\begin{matrix}\frac{d}{d t} x{\left (t \right )}\\\frac{d}{d t} \operatorname{v_{x}}{\left (t \right )}\end{matrix}\right]$$
In [10]:
sol = sympy.solve(model.eqs, list(model.v) + list(model.x.diff(t)))
sol
Out[10]:
$$\left \{ \operatorname{damper.f}{\left (t \right )} : - damper.c \operatorname{v_{x}}{\left (t \right )}, \quad \operatorname{damper.v}{\left (t \right )} : \operatorname{v_{x}}{\left (t \right )}, \quad \operatorname{spring.f}{\left (t \right )} : - spring.k x{\left (t \right )}, \quad \operatorname{spring.x}{\left (t \right )} : x{\left (t \right )}, \quad \frac{d}{d t} \operatorname{v_{x}}{\left (t \right )} : - damper.c \operatorname{v_{x}}{\left (t \right )} - spring.k x{\left (t \right )} - u{\left (t \right )}, \quad \frac{d}{d t} x{\left (t \right )} : \operatorname{v_{x}}{\left (t \right )}\right \}$$
In [11]:
model.x.diff(t).subs(sol)
Out[11]:
$$\left[\begin{matrix}\operatorname{v_{x}}{\left (t \right )}\\- damper.c \operatorname{v_{x}}{\left (t \right )} - spring.k x{\left (t \right )} - u{\left (t \right )}\end{matrix}\right]$$
In [12]:
model.linearize_symbolic()
Out[12]:
$$\left [ \left[\begin{matrix}0 & 1\\- spring.k & - damper.c\end{matrix}\right], \quad \left[\begin{matrix}0\\-1\end{matrix}\right], \quad \left[\begin{matrix}1 & 0\\0 & 1\end{matrix}\right], \quad \left[\begin{matrix}0\\0\end{matrix}\right]\right ]$$
In [13]:
ss = control.ss(*model.linearize())
ss
Out[13]:
A = [[ 0. 1. ]
[-2. -0.2]]
B = [[ 0.]
[-1.]]
C = [[1. 0.]
[0. 1.]]
D = [[0.]
[0.]]
In [14]:
control.ss2tf(ss)
Out[14]:
Input 1 to output 1:
-1
---------------
s^2 + 0.2 s + 2
Input 1 to output 2:
-s
---------------
s^2 + 0.2 s + 2
In [15]:
control.bode(ss[0,0], omega=pl.logspace(-1,1,100));
/home/jgoppert/anaconda3/envs/phd/lib/python3.6/site-packages/control/freqplot.py:153: MatplotlibDeprecationWarning: pyplot.hold is deprecated.
Future behavior will be consistent with the long-time default:
plot commands add elements without first clearing the
Axes and/or Figure.
plt.hold(True);
/home/jgoppert/anaconda3/envs/phd/lib/python3.6/site-packages/matplotlib/__init__.py:805: MatplotlibDeprecationWarning: axes.hold is deprecated. Please remove it from your matplotlibrc and/or style files.
mplDeprecation)
/home/jgoppert/anaconda3/envs/phd/lib/python3.6/site-packages/matplotlib/rcsetup.py:155: MatplotlibDeprecationWarning: axes.hold is deprecated, will be removed in 3.0
mplDeprecation)
/home/jgoppert/anaconda3/envs/phd/lib/python3.6/site-packages/control/freqplot.py:163: MatplotlibDeprecationWarning: pyplot.hold is deprecated.
Future behavior will be consistent with the long-time default:
plot commands add elements without first clearing the
Axes and/or Figure.
plt.hold(True);
In [16]:
control.rlocus(ss[0,0], klist=pl.logspace(-2,1,1000));
pl.grid()
In [17]:
model.f
Out[17]:
$$\left[\begin{matrix}\operatorname{v_{x}}{\left (t \right )}\\- damper.c \operatorname{v_{x}}{\left (t \right )} - spring.k x{\left (t \right )} - u{\left (t \right )}\end{matrix}\right]$$
In [18]:
import sympy
ss_sub = {}
ss_sub.update({model.x[i]: sympy.DeferredVector('x')[i] for i in range(len(model.x))})
#ss_sub.update({model.y[i]: sympy.DeferredVector('y')[i] for i in range(len(model.y))})
ss_sub
Out[18]:
$$\left \{ \operatorname{v_{x}}{\left (t \right )} : x[1], \quad x{\left (t \right )} : x[0]\right \}$$
In [19]:
res = model.simulate(t0=10, tf=20, dt=0.001)
In [20]:
pl.plot(res['t'], res['x'])
pl.grid()
pl.xlabel('t, sec')
pl.ylabel('x')
Out[20]:
Text(0,0.5,'x')
In [21]:
sympy.init_printing()
model.f
Out[21]:
$$\left[\begin{matrix}\operatorname{v_{x}}{\left (t \right )}\\- damper.c \operatorname{v_{x}}{\left (t \right )} - spring.k x{\left (t \right )} - u{\left (t \right )}\end{matrix}\right]$$
In [22]:
model.f.jacobian(model.x)
Out[22]:
$$\left[\begin{matrix}0 & 1\\- spring.k & - damper.c\end{matrix}\right]$$
Using the energy function as a candidate function, we can show that the system has global exponential stability (GES).
In [23]:
v_x, spring_k, t, V, x = sympy.symbols('v_x, spring.k, t, V, x')
V_eq = (v_x(t)**2 + spring_k*x(t)**2)/2
sympy.Eq(V, V_eq)
Out[23]:
$$V = \frac{spring.k}{2} x^{2}{\left (t \right )} + \frac{1}{2} \operatorname{v_{x}}^{2}{\left (t \right )}$$
In [24]:
dV_dt_eq = (sympy.Matrix([V_eq]).jacobian(model.x).dot(model.f)).simplify()
sympy.Eq(V(t).diff(t), dV_dt_eq)
Out[24]:
$$\frac{d}{d t} V{\left (t \right )} = - \left(damper.c \operatorname{v_{x}}{\left (t \right )} + u{\left (t \right )}\right) \operatorname{v_{x}}{\left (t \right )}$$
In [25]:
V_lam = sympy.lambdify(model.x, V_eq.subs(model.p0))
V_vals = [V_lam(*x) for x in res['x']]
pl.plot(res['t'], V_vals)
pl.xlabel('t')
pl.ylabel('V')
pl.grid()
In [26]:
model.f.jacobian(model.p)
Out[26]:
$$\left[\begin{matrix}0 & 0\\- x{\left (t \right )} & - \operatorname{v_{x}}{\left (t \right )}\end{matrix}\right]$$
Content source: jgoppert/pymola
Similar notebooks: