forked from yeus/fmupy
-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathscipy_intergrate.py
65 lines (51 loc) · 1.87 KB
/
scipy_intergrate.py
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
#!/usr/bin/python
# -*- coding: utf-8 -*-
from fmusim import *
import matplotlib.pyplot as plt
import scipy
from scipy.integrate import ode
#myfmu = fmu("./FMU/om/buildingblocks_verosim_basic.fmu", logging=False) #turn logging off for faster calculation
#myfmu = fmu("./FMU/om/satcomponents_blocks_noisetest.fmu", logging=False) #turn logging off for faster calculation
myfmu = fmu("./ibossmo_components_Examples_interfacecomplete_with_dcdc.fmu", logging=False) #turn logging off for faster calculation
contvars=list(myfmu.getContinuousVariables().values())
statenames=list(myfmu.getStateNames().values())
#print(contvars, statenames)
if statenames != []: states = True
else: states = False
varnames = contvars#['imu.y.[1]','imu.y.[2]','imu.y.[3]']
f = myfmu.f #load right-hand-side functio of ode of the FMU
t0,t1 = 0.0,10.0
dt=1.0
y0, status, eventInfo = myfmu.initialize(0.0)
print(myfmu.getStateNames())
res = [[myfmu.getValue(varname) for varname in varnames]]
#file = open("tmp", "w+")
####################################
##start scipy integration
#initialize scipy ode solver
r = ode(f).set_integrator('dopri5') #more methods: http://docs.scipy.org/doc/scipy/reference/generated/scipy.integrate.ode.html#scipy.integrate.ode
r.set_initial_value(y0, t0)
t , y = [t0],[y0]
res = [[myfmu.getValue(varname) for varname in varnames]]
for tn in np.arange(t0,t1+dt,dt):
if states: r.integrate(tn)
y += [r.y]
t += [r.t]
#if not r.successful(): break
myfmu.fmiCompletedIntegratorStep()
res += [[myfmu.getValue(varname) for varname in varnames]]
print(res)
#####################################
##plot result
y = np.array(y)
res = np.array(res)
def plot():
#plot state variables
#for i in y:
#plt.plot(t,y)
#plot only specific variables
for y, var in zip(res.T,varnames):
plt.plot(t,y,label=var.replace('_',''))
plt.legend()
plt.show()
plot()