DyMat を使ってシミュレーション結果を読み込む

DyMat は Dymola や　OpenModelica, JModelica.org などの結果ファイル(matファイル)を読み込んだり、他のファイル形式にエクスポートすることができるライブラリです。これを使って JModelica.org によるシミュレーション結果ファイルのデータをプロットします。DyMat の詳細な使用方法は、Dymat User Manual Dy Mat Guideに記載されています。

まず、mat ファイルのリストを表示します。

In [1]:
%ls *.mat

CSTR_CSTR_Init_Optimization_result.mat
CSTR_CSTR_result.mat
ClassExample1_MassA_result.mat
ClassExample2_MassB_result.mat
EngineV6_analytic_with_input_result.mat
Modelica_Mechanics_Rotational_Examples_First_result.mat
RLC_Circuit_result.mat
VDP_result.mat

Van der Pole oscillator のシミュレーション結果 VDP_resullt.mat を読み込みます。

In [2]:
import DyMat
d = DyMat.DyMatFile('VDP_result.mat')

すべての変数名を表示します。

In [3]:
for name in d.names():
    print(name)

_use_Brent_in_1d
_rescale_after_singular_jac
x2_0
x1_0
_nle_solver_use_last_integrator_step
_block_solver_profiling
_nle_active_bounds_mode
_rescale_each_step
_events_tol_factor
_block_solver_experimental_mode
_nle_solver_default_tol
_nle_solver_check_jac_cond
_use_newton_for_brent
_cs_solver
_block_jacobian_check_tol
_cs_step_size
_events_default_tol
_block_jacobian_check
_iteration_variable_scaling
_nle_jacobian_update_mode
_nle_solver_max_iter
der(x1)
_nle_solver_step_limit_factor
_nle_solver_use_nominals_as_fallback
_nle_solver_min_residual_scaling_factor
_nle_jacobian_calculation_mode
_nle_solver_exit_criterion
_time_events_default_tol
_residual_equation_scaling
_nle_jacobian_finite_difference_delta
_log_level
x2
_nle_brent_ignore_error
x1
_nle_solver_regularization_tolerance
_nle_solver_max_residual_scaling_factor
_nle_solver_min_tol
_use_jacobian_equilibration
_cs_experimental_mode
_runtime_log_to_file
u
_enforce_bounds
_nle_solver_tol_factor
_cs_rel_tol
der(x2)
_nle_solver_max_iter_no_jacobian

アンダーバーから始まる変数名を除いて表示します。

In [4]:
for name in d.names():
    if name.find('_') != 0:
        print(name)

x2_0
x1_0
der(x1)
x2
x1
u
der(x2)

プロットするデータを抽出します。

In [5]:
t = d.abscissa('x1', valuesOnly=True)
x1 = d.data('x1')
x2 = d.data('x2')

データをプロットします。

In [6]:
%matplotlib notebook
import matplotlib.pyplot as plt
plt.figure(1)
plt.plot(t, x1, t, x2)
plt.legend(( 'x1', 'x2'))
plt.title('Van der Pol scillator.')
plt.ylabel('Angle (rad)')
plt.xlabel('Time (s)')
Out[6]:
Text(0.5,0,u'Time (s)')
In [ ]: